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ABSTRACT 



This thesis investigates the problem of automatic 
transcription of the morse signal, and describes and 
documents several approaches to filtering, processing, and 
decoding it for transcription. The baseband signal is first 
modeled as a modified random telegraph wave. A discrete 
Kalman filter and a linear smoother are then used to process 
the demodulated signal in order to gain a measure of the 
effectiveness and applicability of this model. It is shown 
experimentally that this model and processing yield a 
significant reduction in the transcription error rate. Next, 
a Viterbi decoder algorithm based on a simple Markov model 
of the code is programmed and tested. Finally, the base- 
band signal model is incorporated in a more general model 
for pre-detection Kalman filtering. It is shown that this 
filter permits acceptable recovery of morse signals whose 
average signal-to-noise ratio is as low as -14 dB in a 2 kHz 
bandwidth. 
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I. 



INTRODUCTION 



Economic inflation and the national commitment to the 
all-volunteer Armed Forces concept have combined to produce 
unprecedented increases in costs of both manpower and weapon 
systems in recent years. The public's keen awareness of 
these higher costs, together with less than enthusiastic 
support of defense programs in general, has caused Congress 
to be reluctant to authorize increases in defense expenditures. 
Thus it has become necessary to reduce the number of armed 
forces personnel in order to keep defense expenditures within 
authorized limits. This reduction has had the effect of 
intensifying the development of mechanization of appropriate 
manual tasks on a broad front. 

Signal surveillance conducted by the armed forces, 
recognized as an essential and integral part of intelligent 
tactical and strategic planning, is one such area where 
automation is receiving increased attention and support. In 
particular, the human operator has long been relied upon to 
provide the necessary manual transcription of manual morse 
circuits under surveillance. Because of the reduced manpower 
levels, this surveillance and transcription must be tranferred 
to mechanized equipment if this source of intelligence is to 
remain timely and effective. 

This thesis investigates the problem of automatic trans- 
cription of the morse signal, and describes and documents 
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several approaches to filtering, processing, and decoding 
it for transcription. The baseband morse signal process is 
first modeled as a modified random telegraph wave. A discrete 
Kalman filter and a linear smoother are then used to process 
the demodulated signal in order to gain a measure of the 
effectiveness and applicability of this model. It is shown 
experimentally that this model and processing yield a 
significant reduction in the transcription error rate. Next, 
a Viterbi decoder algorithm based on a simple Markov model 
of the code is programmed and tested. Finally, a more general 
model of the signal process, incorporating the baseband model, 
is used to design and implement a pre-detection Kalman filter. 
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II. PROBLEM DESCRIPTION 



A. THE MANUAL MORSE SIGNAL PROCESS 

It is assumed that the reader is familiar with the 
manual morse signal, its pecularities , vagaries, and 
uncertainties, and with current methods of transcription. 

To formalize the discussion, however, certain definitions 
of the terms used and a brief description of the signal are 
in order. 

As used throughout this report, the term morse signal 
refers to International Morse Code, sent manually by key, 
manual "bug", or electronic keyer. The problem of tran- 
scribing keyboard morse , that which is sent automatically 
with standard parameters, will not be considered, although 
certain results are applicable. The baseband morse signal 
is the output of the keyer and is represented by the logic 
levels "0" and "1", corresponding to the states "key up" 
and "key down." The five characters of the international 
morse code are identified as: dot , dash , element-space , 

letter-space , and word-space . The term element ^ refers to 
the standard time unit of the code; its actual duration in 
seconds will of course vary with sending speed. Standard 
morse code consists of the character lengths shown in 
Table I. 
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Sometimes the terms baud and bit are used. 
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TABLE I 



STANDARD MORSE CHARACTERS 



Character 



'Symbol 



Duration (in elements) 



Dot 



1 



Dash 



3 



Element-space 



1 



Letter-space 



3 



Word-space 



w 



7 



The standard word (including word-space) in morse 
communication is 50 elements in length. Thus the element 
duration in seconds for a given sending speed may be 
calculated as 6/5 times the reciprocal of the speed in words 
per minute. The author is unaware of any generally accepted 
standard for the bandwidth of the baseband morse signal; it 
was found to be convenient to express the upper limit of the 
bandwidth as three times the reciprocal of the element 
duration. Thus a code speed of 36 wpm has a bandwidth of 
90 Hz. 

An actual (as opposed to standard) morse signal, as 
those familiar with the problem are aware, may exhibit 
quite a wide variation from standard code in character 
duration, speed variability, and consistency of element 
duration. Since these variations are often unique to the 
particular sending operator, and in many cases may depend 
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on the type of traffic being sent as well, it is difficult, 
if not impossible, to describe a "typical" morse signal. 

This variability is illustrated in Figures 1 and 2, which 
are histograms of the character duration of two different 
amateur ratio operators recorded on the air. As can be seen, 
the distributions are different, with the least variability 
in both cases appearing in the dot and element-space 
durations . 

B . SYSTEM CONSTRAINTS 

The signal processor will obviously play an important 
role in the automated manual morse intercept and transcription 
system. In order to fully appreciate the system constraints 
under which the processor was developed, and the context in 
which it is expected to operate, an outline of the integrated 
system is presented. Referring to the system block diagram. 
Figure 3, its three basic components may be briefly described 
as follows. 

1 . Modulation Sorting Subsystem 

This subsystem scans the band and/or frequencies of 
interest and detects the presence of morse signals in a 
(typically) 2 kHz band. Upon detection of a morse signal, a 
separate digitally-tuned receiver is automatically tuned to 
the signal frequency for reception. 

2 . Morse Processor 

This signal processor is currently the only part 
of the system which is not in existing hardware or software. 
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Its basic function is to minimize error probability in the 
face of noise, interference, and uncertain signal parameters. 

3 . Automatic Transcriber 

This component translates the code into letters of 
the alphabet; there are currently several transcribers 
available which have proven effective at adequate signal- 
to-noise ratios with modest signal parameter variation [1],[2]. 
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Figure 1. Character Duration Histogram 
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III. DESIGN OBJECTIVES 



The ultimate goal of the intercept and transcription 
system is to provide output copy with an error rate no 
greater than that which a "good" human operator can provide. 
Such an operator manually transcribing the morse signal can 
adapt rapidly to changing signal parameters and has little 
trouble distinguishing dots from dashes even if the sender's 
keying is far from perfect. Additionally, he can adapt 
readily to the noise and interference environment and reliably 
copy a signal in the presence of numerous other morse and 
non-morse signals. 

Specific operator performance data were not available to 
the author, thus only broad design objectives were formulated 
based on a limited number of both subjective and experimental 
data obtained using amateur radio operators as subjects. 

Random letter sequences were sent using the Pickering model 
KB-1 morse keyboard to key a signal generator at an audio 
frequency selected by the subject. Noise was added to the 
audio signal and the SNR in each bandwidth used was recorded. 
The results, summarized in Table II, tabulate error rate 
versus SNR in the bandwidth used. Also shown is the SNR in 
the signal bandwidth as previously defined in Section II. A. 

The conclusions drawn from these data is that a good 
operator can copy reasonably well down to -13 dB SNR in a 
2 kHz bandwidth. Although the insertion of a 100 Hz bandpass 
filter raises the SNR to 0 dB , the relative invariance of 
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TABLE II 



OPERATOR PERFORMANCE DATA 



(a) Speed: 35 wpm 



BW 


SNR 


SNR 


ERROR RATE (%) 


COMMENTS 




(dB) 


(dB) 










(in given BW) 


(in signal BW) 


OP1 


OP2 




2 kHz 


-13 


0 


12 


15 


difficult and 


200 Hz 


-3 


0 


11 


13 


fatiguing 


100 Hz 


0 


0 


12 


14 


BW too narrow 


2 kHz 


-10 


3 


6 


10 


fatiguing 


200 Hz 


0 


3 


11 


9 


BW too narrow 


100 Hz 


3 


3 


11 


10 




2 kHz 


-7 


6 


1 


2 


relatively easy 


200 Hz 


3 


6 


0 


1 


prefei wider BW 


100 Hz 


6 


6 


2 


2 





(b) Speed: 25 wpm 



BW 


SNR 


SNR 


ERROR RATE (%) 


COMMENTS 




(dB) 


(dB) 










(in given BW) 


(in signal BW) 


OP1 


OP2 




2 kHz 


-13 


2 


6 


5 


difficult and 


200 Hz 


-3 


2 


7 


6 


fatiguing 


100 Hz 


0 


2 


8 


5 




2 kHz 


-10 


5 


2 


1 


relatively easy 


200 Hz 


0 


5 


1 


1 


but still mildly 


100 Hz 


3 


5 


2 


2 


fatiguing 


2 kHz 


-7 


8 


1 


0 


easy enough 


200 Hz 


3 


8 


0 


1 




100 Hz 


6 


8 


0 


2 
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error rates with changes in filter bandwidth indicates that 
the ear performs the necessary filtering. Strict concentration 
is required, however, at this low SNR, and the test operators 
stated they would not attempt to copy such a signal unless 
strongly motivated. 

Using the previously defined signal bandwidth for each 
speed, the results may be summarized as follows: An operator 

can provide copy with a 10-15% error rate with an SNR of 
approximately 0 dB in the signal bandwidth; at 3 dB , the 
error rate is 5-10%; and at 6 dB, copy is practically perfect. 

Based on these results, it seems reasonable to assert 
that a typical field operator, faced with searching for and 
copying morse traffic eight hours a day, would not be inclined 
to copy signals much below 6 dB SNR in the signal bandwidth 
unless absolutely required. Thus the system designer must 
select from two costly alternatives: 1) If he designs the 

system to perform as well as the good operator is_ able to 
perform, the automated system will reliably receive a large 
percentage of signals encountered on the air, but it is likely 
to be complex and expensive. 2) On the other hand, if he 
designs the system to perform as well as a typical operator 
probably performs, the automated system will be cheaper, but 
the remaining operators who must copy the low SNR signals 
which are not machine transcribable may become too fatigued 
to be efficient, leaving the overall man/machine surveillance 
system less effective than the existing manual system. 

Such design considerations are beyond the scope of this 
thesis; using the "good" operator as a criterion, the 
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ultimate system design objectives may be broadly stated as 
follows : 

1) With an error rate of 10% or less, recover and 
decode morse signals whose SNR in a 2 kHz bandwidth 
is on the order of -10 dB, using standard code with 
additive white gaussian noise and no interference. 

2) Track the time-varying statistics of character 
lengths in order to enable the transcriber to 
translate the code properly. 
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IV. PROCESSOR DESIGN PHILOSOPHY 



The first step in the processor design was to model the 
morse code as a random telegraph wave with non-stationary 
transition probabilities. Using this model, several 
increasingly complex processing methods were implemented, 
and the processing gain of each stage was determined. First 
a Kalman filter was designed to filter the demodulated 
output of a square-law detector with and without narrowband 
analog IF filtering. Next a smoothing algorithm was added 
to determine the effectiveness over Kalman filtering alone. 
Finally a maximum a posteriori (MAP) estimator of the code 
characters, using the smoothed output for likelihood 
calculation, and using the Viterbi algorithm for processing, 
was programmed and tested. 

After determination of the error-reduction effectiveness 
of each of these processing stages, a more general model 
of the signal process was used to design a Kalman filter for 
pre-detection filtering. The objective was to determine 
whether or not such filtering yielded any advantage over the 
simpler demodulation/post-detection filter approach. A block 
diagram of the various stages is shown in Figure 4. Sections 
V through VII present a theoretical basis for each of the 
processing stages, followed by a presentation and discussion 
of experimental results . 
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Figure 4. Processor Design Stages 






V. BASEBAND MODEL AND PROCESSOR 



The baseband morse signal may be modeled as 



x (k+1) = x (k) + w (k) 



(V-l) 



where x = key state (0 or 1) , and w(k) is a random forcing 
function descriptive of the morse keying process. At the 
output of the demodulator, the signal is observed as 



z(k) = x(k) + v(k) 



where z (k) is the observed value and v(k) is the additive 
noise. This model gives rise to the following scalar 
Kalman filter algorithm [3], [4]: 



GOO = v(k l. ^ . ) ... 

V (k | k-1) + R 



(gain) 



V (k | k) = [1 - G(k) ]V(k|k-l) 



(estimation variance) 



V(k+l|k) = V(k|k) + Q (k) 



(prediction variance) 



x (k | k) = x (k | k-1) + G(k)[z(k) - x(k|k-l)] 

(estimation) 



x (k+1 | k) = x (k | k) 



(prediction) 



where : 



A «* . 

x(k) = estimate of x at time k 

R = variance of the observation noise, v(k) 

Q(k) = variance of the random forcing function, w(k). 

Since this algorithm implies knowledge of the variances R 
and Q, procedures for their estimation are required. In 
order to keep the filter algorithm itself as simple as 
possible, estimates of Q and R were made independently, and 
used by the filter algorithm as if these were the true values. 

A. ALGORITHM FOR ESTIMATING THE RANDOM 

FORCING FUNCTION VARIANCE 

The random forcing function w(k) is descriptive of the 
on-off keying process, i.e., it may be thought of as the 
process which causes the transitions in the x state from 
0 to 1 and from 1 to 0 . Referring to the signal model 
equation (V-l) , the keying process has the following 
interpretation: If x(k) = 0 and w(k) = 0, then x(k+l) = 0 

and x remains in the space condition. If, on the other hand, 
w(k) = 1, x shifts (at time k+1) from the space condition to 
the mark condition. This process is illustrated in Figure 5. 
Since the probability of occurrence of a transition is 
dependent on the time duration since the last transition, 
this probability is non-stationary . A proper description of 
the transition probabilities at each time k, then, must 
necessarily be time dependent and conditioned on the element 
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Figure 5 . Illustration of Signal Model Process 
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duration and the present state. The following probabili ties 
are therefore required: 



P (w(k) =0 | k,x(k) =1) = Pr [remain in mark condition] (V-2) 

P(w(k)=l|k,x(k)=l) = 0 ( V- 3) 

P (w(k) =-l |k,x(k) =1) = Pr [transition from mark to space] (V-4) 

P (w(k) =0 |k,x(k) =0) = Pr [remain in space condition] ( V— 5 ) 



P (w(k) =1 | k,x(k) =0) = Pr [transition from space to mark] (V-6) 

P(w(k)=-l|k,x(k)=0) = 0 (V-7) 



Probabilities (V-3) and (V-7) are identically zero since 
state values less than 0 and greater than 1 are not allowable. 
The remaining probabilities are dependent on the distributions 
of dot, dash, element-space, letter-space, and word-space 
durations of the particular morse signal being received, and 
are dependent in a Markov sense on the previous character. 

The Markov nature of the code character transitions was not 
taken advantage of in the filtering process. Additionally 
it was assumed that a particular operator's character 
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durations are all distributed with uniform density, that 
the dot and element space have the same density, and that 
the dash and letter space have the same density. The 
assumed densities, as shown in Figure 6, presume for the 
present that the mean values T, T^, and the parameters £, d 
are either known or have been determined in some manner. 

Although the assumption of uniform densities for 
the character durations of any particular sending operator 
is probably not strictly justifiable, neither is the 
assumption of any other well-known density, such as a 
gaussian or exponential density. The most likely candidate 
for properly modeling these duration distributions may be 
a gaussian-like density with the tails truncated at suitable 
values. The complexity of estimating the parameters of 
such a density for a particular received signal, together 
with the computational burden of evaluating the error 
function, erf(t), for the probability calculations (V-2) 
through (V-7) at each sample point, motivated the selection 
of the uniform density. The resulting probability computations 
are relatively simple and straightforward. 



a) Dot/Element-space Duration Density 




b) Dash/Letter-space Duration Density 




Figure 6. Character Duration Densities 
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1. Mark Transition Probabilities 

A description of probabilities (V-2) and (V-4) may 
be obtained by first conditioning on the dot and dash 
probabilities, and noting that if x(k) is still in the mark 
condition after T+A , then the mark is known to be a dash. 
Given then that a mark, m, is a dot, and given T and A, 
probability (V-2) may be modeled as 



P Q1 (dot) 



A 



T+A 

P(w(k)=0 |k,T,x(k)=l,m=dot) = / Pdot^ < ^ t 



k 



1 



0 < k < T-A 




T-A < k < T+A 



0 



T+A < k 



Similarly, for a dash: 



P Q1 (dash) - P(w(k)=0 |k,T^,x(k)=l,m=dash) 



T d+ d 




1 



0 < k £ T d -d 




T^-d £ k £ T^+d 



0 



T,+d < k 
d — 
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Here k is the time index which counts the number of units 



the signal has been in the mark condition on the present 
mark. The conditioning on m must be removed for 0 _< k £ T+t 
since m is not known during this interval. Thus for k < T+£: 



But the dot and dash probabilities are dependent on the type 
of traffic in the message, i.e. on what language the message 
is in, whether it is plain text or code groups, letters only 
or both letters and numbers, etc. Since the traffic type 
may not be known a priori, equiprobable dots and dashes were 
assumed at this point. Using Pr(dot) = Pr(dash) = 1/2, then, 
equations (V-8) and (V-9) reduce to 



P (w (k) =0 |k,T,x(k)=l) = P 01 (dot)Pr(dot) + P Q1 (dash) Pr (dash) 



(V-8) 



and for k > T+i. 



P(w(k)=0 |k,T d ,x(k)=l) = P 01 (dash) 



(V-9) 



1 



0 < k < T-S, 



T-k 3 



4 z 4 



T-Z < k < T+5, 



P(w(k)=o|k,T,T d ,x(k)=l) 



1 



T+& < k £ T (j“ d 




T d ~d <_ k <_ T d +d 



0 



T d +d < k 
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Probability (V-6) follows immediately as 



P(w(k)=-l|k,T,T d ,x(k)=l) = 1 - P(w(k)=0 |k,T,T d ,x(k)=l) . 

The expression for probability (V-2) as a function of time, 
k, conditioned on the present state and dot and dash 
distributions is sketched in Figure 7. 

2 . Space Transition Probabilities 

An appropriate description for the space duration 
probabilities (V-4) and (V-7) is derived similarly, first by 
conditioning on a particular space and then removing the 
appropriate conditioning by using the relative frequencies 
of each space. Given that a particular space, s, is an 
element-space, then. 



1 0 < j <■ T-£ 

P (w( j) =0 | x( j ) =0 ,s=elem,T) = + J T-£ £ j £ T+£ 

0 T+£ < j 

(V-10) 

where j is the time index which counts the number of units 
the signal has been in the space condition on the present 
space; and the element length, given T, is assumed to have 
the same uniform density as the dot length. 
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Probability 
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Similarly for a letter-space: 



P(w( j)=0 | x(k) =0 ,s=ltr-sp,T ,) = 



o < j < T d -d 



i 

2d + 2 



T d " d 1 3 1 T d +d 



(V-ll) 



T d + d < j 



where j is the same index as above and the letter-space 
length, given T^, has the same density as the dash- length. 

For the word-space it was decided to use an 
exponential model, since after about 5T units, the word-space 
is (in actual practice) about equally likely to end at any 
time: 

1 0 < j < 5T 



P (w ( j ) =0 | x( j ) =0 , s=word-sp ,T) = 



-( 



e 



j~5T . 
2T 1 



5T < j < 00 



(V-12 ) 
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The conditioning on s may be removed by observing that if 
j > T^+d, the space must be a word-space; if j > T+£, then 
the space is a word-space or letter-space. In equiprobable 
letter traffic, the needed probabilities are: 

Pr (s=elem-space) =12/17 
Pr (s=letter-space) = 4/17 
Pr (s=word-space) = 1/17 

Applying these probabilities then to (V-10) , (V-ll) , and 
(V-12) yields the desired expression: 



P(w( j)=0 | j,T,T d ,x(j)=0) 



1 



0 < j < T-& 




T-£ < j < T+ Jt 



1 



T+Jl < j <_ T d ~d 




T,-d < j < T ,+d 
d J — d 



1 



T ,+d < j < 5T 
d J — 




5T < j 



e 
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3. Calculation of Variance 



The variance of w(k) is now easily calculated at 
each time, k, as follows: 

Let 



P k (0,l) = Pr (w (k) =0 | x (k) =1) 

P k (-l,l) = Pr (w (k) =-l | x (k) =1) 

Pj(0,0) = Pr (w( j)=0 | x ( j ) =0 ) 

Pj(l,0) = Pr (w( j ) =1 | x ( j ) =0) . 

Then, if x(k) = 1, 

Q(k) = (0) 2 -P k (0,l) + (-l) 2 -P k (-l,l) - [E (w (k) ) ] 2 
= P k (-l,l) - t-P k (-l,l) ] 2 
= P k (-l,l) .P k (0,l) . 

Similarly, if x(j) = 0, 

Q ( j) = Pj(0,0) ^(1,0) . 

Since, at the receiver, the true state x(k) is not 
known, correct estimates of Q are dependent on correct 
estimates of the state, x(k) . Thus, it is expected that at 
some SNR, incorrect estimates of x(k) will cause the 
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estimates of Q to be erroneous enough to force incorrect 
estimates of x(k+l) , causing a runaway condition to develop 
and yield the receiver worthless. This SNR at which runaway 
develops was determined experimentally, and found not to be 
a serious problem. 

B. CHARACTER DISTRIBUTION ESTIMATION 

It was assumed that the character distributions are 
appropriately described by uniform densities with known 
parameters (expect for mean values) . More sophisticated 
methods of distribution and/or parameter estimation were 
discarded in favor of simplicity. Experience suggests that 
even a sloppy sender will usually not exceed an £/T or d/T^ 
ratio of about 1/3. Thus, once estimates of T and T^ are 
obtained, i and d can be determined from this assumption 
unless these ratios are known in advance. The mean values 
T and T^ of dot or element-space and dash or letter-space, 
respectively, were estimated by measuring the character 
durations and sequentially averaging the appropriate duration. 
Thus T is the mean value of the dot and element-space durations 
and T^ is the mean value of the dash and letter-space durations 
The threshold for deciding which set of measurements a 
particular length belongs to was set at the halfway-point 
between dot and dash lengths. The algorithm is as follows: 
Initially specify: 

t^ = shortest dot-duration expected, 

t 2 = longest dot-duration expected, 
with t 2 £ 3t^. 
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1 . 



(a) Initialize estimate of T as (t 1 +t 2 )/2. 

(b) Initialize estimate of as (3t^+3t 2 )/2. 

2. Measure the dura tion J ~of each mark and space. 

3. If the measured value is less than (3t^+t 2 )/2, then 
identify it as a dot/element duration, T^. 

4. If the measured value is larger than (3t 1 +t 2 )/2, then 
identify it as a dash/letter-space duration, T 2 . 

5. Estimate the means recursively by 



(a) 


T(k) 


= T(k-l) + i [T x - T(k-l) ] . 


(b) 


V k > 


‘ v*- 1 * + E [T 2 - V k -in 



C. OBSERVATION NOISE VARIANCE 

Although the noise power used in each of the tests was 
known and could have been specified initially, it was 
decided to estimate this parameter in order to better 
simulate an operational environment where the noise power 
is not known a priori. Although the method is purely 
intuitive and without a valid theoretical basis, reasonably 
good results were obtained. 

At the output of the square-law demodulator, the noise 
is no longer gaussian and is correlated with the signal. 
Proceeding, however, as if the signal and noise were not 
correlated, the noise variance R can be obtained by sub- 
tracting the morse signal power from the total received 
(demodulated) signal power: 
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R = 



A 





A 

where V is obtained recursively by the following: 
V k) = v z (k-l) + i [(z-i z ) 2 - V z (k-1) 3 



and 



y = mean value of the received (demodulated) signal 
z 

A 

a = estimate of demodulated signal amplitude. 



The parameters y and a were also estimated on line, although 

z 

they too are known a priori for test purposes. The estimator 

A 

algorithm for a is 



1) a x (k) = a 1 (k-l) + £ [z(k) - a^k-1)] 



if z (k) > y z 



(2) a 2 (k) = a 2 (k-1) + ^ [z(k) - a 2 (k-1) ] if z(k) < y. 



(3) a (k) = a^ (k) - a 2 (k) . 



The mean value, y , is simply 



y„ (k) = y (k-1) + £ [z (k) - y (k-1) ] 
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D. SMOOTHING ALGORITHM 



The Kalman-filtered output represents the optimum 
linearly-filtered estimate of the signal amplitude in the 
minimum mean-square-error sense based on the assumed model 
and the efficiency and consistency of the estimated 
statistics. This estimate may be improved, however, by 
smoothing the data, which implies consideration of future 
inputs as well as past inputs [4]. 

The equations for the discrete optimal linear fixed- 
interval smoother, for the scalar case, expressed as a 
combination of a forward- running filter and backward-running 
filter are given by [5], [6] as: 

Smoothed Estimate: 



x (k | N) = P (k |N). [x(k | k) /P f (k |k) + x^ (k | k+1) /P^ (k | k+1) ] 



Estimation Variance: 



P (k |N) 



1 

1/P f (k | k) T 1/P b (k|k+1) ' 



where 



x(k|k) 



P f (k | k) 

x b (k|k+l) 

P b (k|k+1) 



filtered state estimate 
error variance for the forward filter 
(predicted) estimate for the backward filter 
(predicted) error variance, backward filter. 
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These expressions are not computationally suitable, however, 
since P b (k|k+1) is not finite when k = N. Alternate 
expressions [6] which are more amenable to computation are 
given as follows: 



Let 



W(k) 



P f (k|k) 

1 + P f (k|k)/P b (k|k+l) • 



Then 

P(k|N) = [1 -W(k)/P b (k|k+1) ] 2 -P f (k | k) • [W 2 (k)/P b (k|k+1) ] 



and 



1/P^(k|k+1) = 0 when k = N . 



Let W b (k|k+1) = Xj^ (k | k+1) /P fe (k | k+1) . 



Then 



x (k |n) 



x (k | k) 

[1 + P f (k|k)/P b (k|k+l) ] + P(k|N)W b (k|k+l) 



and 



P (k |N) W b (k | k+1) = 0 when k = N . 
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These equations then represent the smoothing algorithm. 

As noted above, forward-filtered estimates must be stored 
until the same data can be backward-filtered and combined 
to produce the smoothed state estimate. 

The smoothed state estimate then is the best (linear) 
estimate possible for the given model and parameter estimators, 
and any improvement in error rate must be derived from giving 
more probabilistic structure to the signal model. 

E. IMPLEMENTATION OF FILTER AND SMOOTHER 

A diagram of the test signal generation method is shown 
in Figure 8. The square-law demodulator was selected simply 
for ease of implementation since an analog squarer and 
appropriate filters were readily available and easily inter- 
faced with the analog-to-digital converter. The 100 Hz 
low-pass filter permits recovery of morse signals of 
approximately 35 wpm or less. The signal-to-noise ratios 
used throughout this report are average (pre-detection) 
signal-to-noise power ratios and not pulse signal power to 
noise power ratios. More specifically, the average signal 
power in the morse signal is defined as 

2T 

p s = 2^r 0 f Cs(t) ]2 dt 

where the interval [0,2T] is a dot and element-space. This 
expression reduces to 
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T 



2 dt + ,p= f 0 dt 

T 



P 



s 




which is 3 dB less than the pulse signal power. 

The noise signal is taken from a calibrated white gaussian 
noise source of 1 volt rms. 

The filter and smoother, along with the auxiliary 
estimation algorithms, were coded in Fortran and implemented 
on the XDS-9300 computer interfaced with the CI-5000 analog 
computer for analog filtering and for D/A and A/D conversion. 
The sampling rate was 500 samples per second, and the value 
of N for smoothing was chosen to be 250 samples. The sampled 
test signal was recorded on tape for subsequent processing, 
since the processing required approximately 4 seconds for 
1 second of data. The test signal runs consisted of the 
following : 

1. Perfect code AR sequence and random letter sequence at 
speeds of 35, 30, and 25 wpm each with a signal- to-noise 
ratio of 6, 5, 4, 3, 2, and 1 dB in a 2 kHz BW with no 
pre-detection analog filtering. 

2. Perfect code AR sequence and random letter sequence at 
speeds of 35, 30, and 25 wpm each with a signal-to-noise 
ratio of -7, -8, -9, -10, -11, -12 dB in a 2 kHz BW with 
a 100 Hz pre-detection analog bandpass filter. 
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3. Sloppy code AR sequence at nominal speeds of 25 and 30 
wpm each with a signal-to-noise ratio of 6, 5, 4, 3, 2, 



1 dB in a 2 kHz BW with no pre-detection filtering. 

4. Sloppy code AR sequence at nominal speeds of 25 and 30 
wpm each with a signal-to-noise ratio of -7, -8, -9, 

-10, -11, -12 dB in a 2 kHz BW with a 100 Hz pre-detection 
analog bandpass filter. 

The output of the processor was recorded on an 8-channel 



strip-chart. 


using a utility D/A conversion routine. The 


outputs were 


as follows: 


Channel 1: 


Recorded input signal 


Channel 2: 


Unprocessed output (input signal thresholded 
at its mean value) 


Channel 3: 


Kalman filtered output 


Channel 4: 


Filtered output thresholded at 0 V. 


Channel 5: 


Smoothed output 


Channel 6: 


Smoothed output thresholded at 0 V. 


Channel 7: 


Option of filter gain or noise variance estimate 


Channel 8 : 


Dot/element-space duration estimate. 



The following figures (9-19) are typical output records, 
showing examples of the test runs for each signal and type 
of sequence. Channels 1-4 are shown in Figure (a) in each 
case with the corresponding channels 5 and 6 shown in Figure 
(b) . An example of the output of channel 7 for the gain 
option is shown in Figure 20, along with the corresponding 
outputs from channels 1, 3, and 4. Figure 21 shows an 
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example of the outputs for the noise variance estimate and 
the dot/element-space duration estimate at a point where the 
input signal changes from 25 to 30 wpm and the signal-to- 
noise ratio changes from 1 dB to 6 dB. In all cases the 
chart speed was 10 mm/sec except in Figure 21 where 5 mm/sec 
was used in order to show the estimates more clearly. The scales 
for channels 1, 2, 4 and 6, and for channel 7 gain option, 
are 5 v./div., with 100 volts corresponding to a variable 

value of 1.0; the scales for channels 3 and 5 are 2 v./div. 

For the channel 7 variance option the scale is 200 mv./div. , 
snd for channel 8 the scale was calibrated at 4 msec/div. 

F. RESULTS OF TESTS 

The outputs of the processor were decoded by hand to 
determine the error rates. Using the estimate of T, a mark 
was decoded as a dot if its duration was 2T or less and as 
a dash otherwise. Similarly, a space was decoded as an 
element-space if its duration was 2T or less, as a letter- 
space if the duration was between 2T and 4T, and as a word 
space otherwise. 

The following data were obtained for each run: 

1. Letter error rate and bit error rate with no processing. 

2. Same error rates with filtering only. 

3. Same error rates with filtering and smoothing. 

A letter error occurs when any transmitted letter is not 
correctly decoded. Only one error per transmitted letter 
is counted; for example, if "A" is decoded as "ET" , one 
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letter error has occurred. On the other hand, if "ET" is 
decoded as "A", two errors have occurred, since neither "E" 
nor "T" was decoded correctly. Word spaces are counted as 
one letter per 7 elements. The sample size in each case was 
approximately 200 letters. A bit error is defined as at 
least one mark-space error occurring within a transmitted 
element duration. Again, only one error per element is counted, 
and the sample size was approximately 200 bits. 

The results of this analysis are presented in Tables III 
through VI. Column 1 lists the signal-to-noise ratio used 
without a 100 Hz pre-detection filter; column 2 lists the 
signal-to-noise ratio used with the 100 Hz filter in place. 

Both the bit and letter error rates for filtering and 
smoothing are tabulated, with the error rates for the 
unprocessed output shown for comparison. In each case, the 
error rates are shown for the typical (random letter) 
sequence and the AR sequence, except in Table VI, where the 
results are for the AR sequence only since no random letter 
sequence for this case was recorded. Table VII shows 
typical hand-translated sequences for each processing stage. 

These results indicate that the Kalman filter and linear 
smoother provide a significant decrease in both bit and 
letter error rates. By using a 100 Hz bandpass pre-detection 
filter, such processing provides a tolerable 10% letter error 
rate on a -7 dB SNR signal as opposed to an unacceptable 
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FIGURE 9a. Filter Output for 35 wpm Sequence, 6 dB SNR 
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FIGURE 9b. Smoothed Output for 35 wpm Sequence, 6 dB SNR 
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FIGURE 10a. Filter Output for 35 wpm Sequence, 5 dB SNR 
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FIGURE 10b. Smoothed Output for 35 wpm Sequence, 5 dB SNR 



Channel 







i — i cm ro ^ 



52 



FIGURE 11a. Filter Output for 35 wpm AR Sequence, 4 dB SNR 
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FIGURE lib. Smoothed Output for 35 wpm AR Sequence, 4 dB SNR 
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FICURE 12a. Filter Output for 35 wpm Typical Sequence, 6 dB SNR 
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FIGURE 13a. Filter Output for 35 wpm Typical Sequence, 4 dB SNR 
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FIGURE 13b. Smoothed Output for 35 wpm Typical Sequence, 4 dB SNR 
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FIGURE 14a. Filter Output for 30 wpm AR Sequence, 5 dB SNR 
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FIGURE 15a. Filter Output for 25 wpm AR Sequence, 5 dB SNR 
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FIGURE 16a. Filter Output for 25 wpin AR Sequence, 3 dB SNR 
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FIGURE 16b. Smoothed Output for 25 wpm AR Sequence, 3 dB SNR 
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FIGURE 17a. Filter Output for 25 wpm Sloppy Code AR Sequence, 6 dB SNR 
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FIGURE 17b. Smoothed Output for 25 wpm Sloppy Code AR Sequence, 6 dB SNR 
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FIGURE 10a. Filter Output for 25 wpm Sloppy Code AR Sequence, 3 dB SNR 
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FIGURE 18b. Smoothed Output for 25 wpm Sloppy Code AR Sequence, 3 dB SNR 
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FIGURE 19a. Filter Output for 25 wpm Sloppy Code AR Sequence, 1 dB SNR 
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FIGURE 19b. Smoothed Output for 25 wpm Sloppy Code AR Sequence, 1 dB SNR 
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FIGURE 20. Filter Gain (Channel 7) 4 dB SNR 
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FIGURE 21. Noise Variance and Estimate of 



ERROR RATES FOR CODE SPEED OF 35 wpm 
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ERROR RATES FOR CODE SPEED OF 30 wpm 
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ERROR RATES FOR CODE SPEED OF 25 wpm 
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ERROR RATES FOR SLOPPY CODE 
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TABLE VII 



TYPICAL TRANSLATED SEQUENCES 



SENT 

UNPROCESSED 

FILTERED 

SMOOTHED 



QA Z WS X E DC R F V T GB YHN UJM 
?? ? E WH XEIED? E FFV T GG E YHN FJM 
?? ? E WH X IEDC E RFV T GG YHN UJM 
?? ? WH X E DC E RFV T GB YHN UJM 



SENT 


QAZ 


WS X E DC 


RFV T 


GB YHN UJM 


UNPROCESSED 


IEQAZ 


PSEX I DCEFF7EA 


GBEYHNEU?? 


FILTERED 


IEQAZ 


?? X I DC 


FFVIT 


GB YHN UJM 


SMOOTHED 


EEQAZ 


?? X T DC 


RF? T 


GB YHN UJM 


SENT 


QAZ 


WSX E E DC 


RFV T 


GB YHN UJM 


UNPROCESSED 


QUAD 


??XEEI I DC E 


LPVET 


PBEYHNE? JME 


FILTERED 


QMZ 


WHXEEI I DY E 


RFV T 


GB YUN F?ME 


SMOOTHED 


QAZ 


WHX I I DC 


RFV T 


GB QUNEUJM 


SENT 


QAZ 


WS X E DC 


RFV T 


GB YHN U 


UNPROCESSED 


Q?EZ E 


LS X EIH ?? 


LFVIA 


GB E YH7EV 


FILTERED 


QAZ E 


WS X EEE DC E 


RFVET 


GB YHN U 


SMOOTHED 


QAZ 


WS X E D? E 


RFV T 


GB YHN U 


SUMMARY: 

total letters sent 


125 


Error Rate 


unprocessed errors 


58 




46% 


filtered 


errors 


31 




25% 


smoothed 


errors 


19 




15% 


(question mark indicates untranslatable 


sequence) 


BW = 2 kHz 


SNR = 5 dB speed = 30 wpm 


perfect code 
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32-34% error rate without processing. Smoothing contributes 
considerably to a reduced error rate down to about -9 dB 
where its effectiveness over filtering alone begins to fall 
off. For the high speed 35 wpm signal, runaway of Q estimation 
occurs at approximately -9 dB , while for the lower speeds 
runaway never really develops until the SNR reaches -12 dB, 
since filtering always provides an improvement in error rate 
even at these low SNR's. 

It was noted that a majority of the errors in the filtered 
and smoothed output result from insertions of isolated dots 
in the letter-space and particularly the word-space separations. 
A possible remedy to this situation is to incorporate the 
Markov structure of the code in the estimation algorithm for 
Q, although an increased susceptibility to runaway may limit 
its effectiveness. 

Since the bit error rate must be on the order of 1% or 
less in order to yield a tolerable 10% letter error rate, 
it was felt that bit error rates of up to approximately 3% 
could possibly be reduced to an acceptable level by use of 
soft-decision Viterbi decoding following the smoothed output. 

If such a reduction were possible, then acceptable error 
rates could be obtained for SNR's down to -9 dB . 

Additionally, pre-detection Kalman filtering and parameter 
estimation would yield a theoretical gain of 3 dB or more 
over non-coherent demodulation. With these improvements in 
the processor, then, the output letter error rate of 10% 
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could be achieved on signals whose SNR is on the order of 
-12 dB in a 2 kHz bandwidth and the originally stated design 
objective would be met. 

Both of these possible improvements were implemented 
separately to determine the effectiveness of each one. The 
Viterbi decoder algorithm used the smoothed output of the 
post-detection Kalman filter as input. A separate program 
was written to implement a pre-detection Kalman filter. 
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VI. VITERBI DECODER 



The Viterbi algorithm [7], [8], as originally formulated, 
is a maximum likelihood (ML) sequence estimation decoding 
algorithm for convolutional codes. It has found application 
in other areas [9], however, and its use has been extended 
to maximum a posteriori (MAP) estimation. It is the MAP 
estimation use which is of importance here. 

A. MAP ESTIMATION 

In order to transform the smoothed output of the 
processor into characters of the morse code, certain 
decision criteria must be formulated. The easiest and most 
obvious way to accomplish decoding is to threshold the 
smoothed output at its mean value and determine the identity 
of marks and spaces on the basis of the measured duration 
of each received character. Such a scheme, however, fails 
to utilize two sources of information which are inherent 
in the smoothed output: 1) Thresholding discards all 

information present in the actual smoothed amplitude estimate, 
and 2) the Markov nature of the character transitions is not 
utilized. 

The decoding problem, then, is to take advantage of this 
additional information to determine the most probable 
sequence of morse characters. The thresholded output may be 
used to make tentative decisions to obtain a specific 
sequence of character outputs z ,...,z , where a particular 
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z^ is either a dot, dash, element-space, or letter-space. 

It is the objective, then, to find a specific sequence of 

transmitted characters x ,...,x which maximizes the 

o n 

probability that x , ...,x n was transmitted, given that 
z , ...,z n was received. Formally, it is desired to find 
the sequence Xj = x Q ,...,x n which gives 

P I ) = maximum 



where z^ is the specific received sequence, and all possible 
sequences of transmitted characters, x, are used to determine 



which specific sequence yields the maximum value. The 
sequence x^ which yields the maximum is then the maximum 
likelihood estimate of the transmitted sequence if it is 
assumed that the transmitted characters are equiprobable; 
it is the maximum a posteriori (or minimum error) estimate 
if the actual probability of transmission of each character 
is utilized [10] . 

In general, it would be necessary to compute and compare 
the probability, p (2£ = 2i£ l^ = £j) / for all possible sequences x^ 
However, if it is assumed that the morse code is a Markov 
source, then the problem of finding the MAP estimate reduces 
to a problem of maximizing a sum and the Viterbi algorithm 
may be used. 

In the following development, a shorthand probability 
notation is used to facilitate writing of probability 
statements. The statement P(x^.|z^) is used to mean 
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P (x^=a^ | z^,=aj ) where the a^,a^ are the characters of the 
code, i.e., dot, dash, element-space, letter-space. No 
confusion should result, since in all cases the notation 
P (u) is intended to mean the probability that u is equal 
to some specific value. 

B. SOURCE MODEL 

A third-order Markov model of the morse code exhibits 
a good deal more probabilistic structure than a first-order 
model, as can be seen by comparing the transition probabilities 
shown in Table VIII. In the interest of simplicity, however, 
it was decided to use a first-order model. The assumption of 
a first order model then means that 



p(x kl x o 



' p * x k' x k-l > 



where x^ is the kth character of a transmitted sequence. 

The transition matrix lists the following transition 
probabi lities : 



P(x k |x k _ 1 ) for each x k = a i 



where 






a l = * 


(dot) 


(Vl-la) 


a 2 = — 


(dash) 


(Vl-lb) 


a 3 = ^ 


(element- space) 


(VI-lc) 


a 4 = ~ 


( letter-space) 


( Vl-ld) 
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TABLE VIII 



MARKOV TRANSITION PROBABILITIES 
FIRST-ORDER MODEL 





• 


/N 






• 


0 


.682 


.318 


0 


A 


.554 


0 


0 


. 446 


- 


.5 


0 


0 


.5 


— 


0 


.684 


.316 


0 



THIRD-ORDER MODEL 





# 


A 


- 


— 


» A • 


0 


.438 


.562 


0 


• A ■*"' 


0 


.615 


.385 


0 


• ^ « 


0 


.923 


.077 


0 




1 ' ■ 1 

0 


.923 


.077 


0 


A • A 


1 .556 

l 


0 


0 


.444 


A • ~ 


.50 


0 


0 


. 50 


A“A 


.538 


0 


0 


.462 


/A — ~ 


.50 


0 


0 


. 50 


a » A 


.50 


0 


0 


.50 


- • - 


.50 


0 


0 


.50 


a“A 


. 583 


0 


0 


.417 


- 


.50 


0 


0 


. 50 


"~A • 


0 


. 571 


.429 


0 


“ “A“ 


0 


. 545 


.455 


0 


— ~ . 


0 


. 923 


.077 


0 




0 


.923 


.077 


0 
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(Word spaces were not considered as separate characters, 
but as combinations of letter-spaces and element-spaces.) 



C. SEQUENCE PROBABILITIES 

Based on the first order model assumption, the probability 
of any particular transmitted sequence of length n is given by 



n-1 

P ( X , x, , . . . , X - ]T P ( x • | x . _ i ) 

u i=0 



(VI-2) 



Then, given an input sequence to the decoder (z , ...,z ^), 

use of Bayes' rule expresses the desired conditional 
probability as 



P (x c , . . . /X n _ 1 | Z Q , . . . , 2 n-1 ) 



P(z 0 ,...,z n _ 1 |x o x n _ 1 ) -P(x o , 



,x , ) 
' n-1 



p(z o z n-l> 



(VI-3) 



which is the probability to be maximized. 

Assuming that the thresholded output is memory less, the 
conditional output sequence probability becomes: 



P ( z 



o' 



z n-l ' x o ' 



x , ) 
n-1 



n-1 

n p(z. | x. ) 
i=0 



(VI-4) 



Although the thresholded output is by no means memory less, 
due to the decision directed nature of the Q estimation 
algorithm in the preceding filter, this assumption is 
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nevertheless made in order to render the computation of the 
necessary probabilities tractable. The dependence of the 
thresholded character output decision on previous and future 
decisions will be removed to a certain extent in Section D, 
justifying the assumption of memory lessness at this point. 

Then, using (VI-2) and (VI-4), and realizing that 
P(z Q ,...,z ^) is constant, maximization of (VI-3) is 

equivalent to maximizing the expression: 

n-1 n-1 

n P(x. |x. ,) n P ( z . | x . ) 
i=0 i=0 

Maximization of this sequence is equivalent to minimizing 
the negative logarithm, since In P(*) is a monotonic function 
of P ( . ) . Thus 

A n ~ 1 

L(x o , . . . ,x n _ 1 ) = - [In P(x i |x i _ 1 ) + In P(z i |x i )] 

is the function to be minimized by the Viterbi algorithm. 

An outline of how the Viterbi algorithm performs this 
minimization is presented in Appendix B. 

D. LIKELIHOOD COMPUTATION 

The likelihoods P(z^|x^) may be computed from the smoothed 
signal amplitude and received character duration as follows. 

Define the following figures of merit for amplitude and 
duration : 
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a^ = amplitude figure of merit for a dot-length duration 
t = time figure of merit for a dot-length duration 
A^ = amplitude figure of merit for a dash- length duration 
= time figure of merit for a dash-length duration 
If these figures of merit are scaled such that their 
values are between 0 and 1, then they may be interpreted 
as follows: 

a m = probability that a mark occurred during a dot-length 
interval. 

t = probability that a mark is a dot, or 

= probability that a space is an element-space. 

A = probability that a mark occurred during a dash-length 
interval . 

= probability that a mark is a dash, or 
= probability that a space is a letter-space. 
Likelihoods, then, may be computed by utilizing these values 
as probabilities. For example, the probability that the 
thresholded output is a dot, given that the input is an 
element-space is simply: 

P (z=a 1 | x=a 3 ) 

= Pr (z=dot | z=mark) *Pr ( z^dash | z=mark) •Pr(z=mark) 
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The likelihoods for each case were determined as above 



and are given as follows : 



P(z=a 1 |x=a.) = a t (1-T ) 

1 1 1 rn m ' nr 

P(z=a~|x=a.) = a (1-t (1-T)) 
2 1 x m m m 



P(z=a 3 |x=a i ) = (l-a m )t in (l-T m ) 



P< Z -a 4 |x=a.> = d-a m )(l - t m (l-T m )) 



for i = 1,3 where the a^ are given by (VI-1); and 



P( Z =a 1 |x=a i ) = ajl - ^(1-t^) 

P(z=a 2 |x=a i ) “ A rn T m < 1 - t m ) 

P( Z =a 3 |x=a i ) = (l-a m )(l- TJl-tJ) 
P(z-a 4 | X =a.) = U-V T n. (1 " t m ) 



for i = 2,4. 

The figures of merit a , A , t , T were obtained by 
y m' m m' m 2 

using the thresholded, smoothed output to make "tentative" 
decisions, and then computing the merit of these decisions. 
The tentative decisions were: 
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1) 

2) 

3) 

4) 



If x 
If x 
If x 
If x 



out 

out 

out 

out 



1 and measured duration <_ 2T, 
1 and measured duration > 2T, 
0 and measured duration _< 2T, 
0 and measured duration > 2T, 



then x k 
then x^ 
then x^, 
then x^. 



dot 



dash 



element- 

space 

letter- 

space 



For Xj, = dot or element-space, the previous T values of 
smoothed output were averaged and the amplitude figure of 
merit was taken as 



0 




1 x g (T) 



a = x (T)/a + 0.5 
m s 



li X S (T) 



> 



a 



2 



1 



A 

x s (T) > | 



where 



X S (T) 



A 

a 



the average of the smoothed output, x g , over the 
previous T seconds. 

amplitude of smoothed output signal. 



The time figures of merit were 



1 



0 < t d < T 



t = 2 - (t , + T) /2T 

m d ' 



T < t , < 3T 
— d — 



0 



3T i 
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0 



o i t a < T d /3 



T m = 1 ' (T rt -t fl )/(0.67T fl ) T d /3 < < I 



d d 



d-d 



T d 1 fc d 



where = measured duration of character duration. 

These functions are sketched in Figure 22. For Xj, = dash 
or letter-space, the previous T d values of x g were averaged 
to obtain the value for A : 



0 



A = x (Tj/a + 0.5 
m s d 




i x s (T a> 



/S 




/V 



a 

2 



1 



/s 




These likelihood computations allow decisions to be made 
by the Viterbi algorithm to determine the most likely 
character on a character-versus-character basis. It was 
necessary, however, to extend these computations to cover 
more complicated situations such as that depicted in Figure 23. 
The thresholded output shown decodes as , although the 

sequence . was actually transmitted. The Viterbi algorithm 
implemented using the above likelihood computations decoded 
the sequence as . , however, since no provision was made 

to account for the non-memory less nature of the thresholded 
output . 
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a) Amplitude Figure of Merit 




b) Duration Figures of Merit 




Figure 22. Sketch of Figures of Merit 
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Smoothed 




Time 



Thresholded 

Output 



1 - 



nme 



Figure 23. Illustration of Signal Requiring 
Modified Likelihood Computation 
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Two modifications were made, therefore, in order to 
account for situations such as this: First, when a dash or 

letter-space was decoded by the thresholded output, this 
character was divided into three equal segments and dot and 
element-space likelihoods were computed for each segment. 
Secondly, when the thresholded output consecutively decoded 
three short characters, then the total duration of all three 
characters was obtained and dash and letter-space likelihoods 
were computed for the total interval. 

E. IMPLEMENTATION AND RESULTS 

The algorithm was coded in Fortran and run as a 
subroutine of the original Kalman filter and smoothing 
programs. Sequences of 9-11 elements in duration were 
decoded. In order to determine the contribution of context 
information supplied by the first order transition 
probabilities, the algorithm was also run with equiprobable 
transition probabilities, i.e. as a ML sequence estimator. 

The results are summarized in Tables IX and X. Again 
both bit and letter error rates are presented with the hand- 
decoded smoothed output repeated here for ease of comparison. 
Table XI shows a comparison of bit error rates for ML and 
MAP estimation showing a slight improvement provided by the 
MAP estimation, taking advantage of transition probabilities. 

The reduction in error rate, although not as great as 
hoped for, shows that the thresholded output bit error rate 
can be improved significantly at the higher SNR's, but 
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TABLE IX 



VITERBI DECODER ERROR RATES - 35 wpm 



AR SEQUENCE, PERFECT CODE: 



SNR (2 kHz) 


SMOOTHED OUTPUT 


VITERBI 


OUTPUT 


IMPROVEMENT RATIO 


NO PKE-FILIER 


BIT 


LTR 


BIT 


LTR 


IN BIT ERROR RATE 


(db) 


(%) 


(%) 


(%) 


(%) 




6 


0.83' 


9 


.35 


3 


2.4 


5 


3.0 


27 


2.5 


23 


1.2 


4 


6.7 


49 


6.8 


49 


- 


3 


35 


96 


35 


100 


_ 



TYPICAL SEQUENCE, 


PERFECT 


CODE : 






SNR (2 kHz) 


SMOOTHED OUTPUT 


VITERBI 


OUTPUT 


IMPROVEMENT RATIO 
IN BIT ERROR RATE 


NO PRE-FILTER 


BIT 


LTR 


BIT 


LTR 




(db) 


(%) 


(%) 


(%) 


(%) 




6 


0.70 


8.0 


0.33 


3 


2.1 


5 


1.2 


15 


0.90 


11 


1.3 


4 


4.5 


34 


4.4 


34 


- 


3 


30 


92 


30 


96 


_ 
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TABLE X 

VITERBI DECODER ERROR RATES - 25 wpm 
AR SEQUENCE, PERFECT CODE: 



SNR (2kHz) 


SMOOTHED OUTPUT 


VITERBI 


OUTPUT 


IMPROVEMENT RATIO 


NO PRE-FILTER 


BIT 


LTR 


BIT 


LTR 


IN BIT ERROR RATE 


(db) 


(%) 


(%) 


(%) 


(%) 




6 


.77 


10 


0.31 


2 


2.5 


5 


1.3 


17 


1.0 


13 


1.3 


4 


1.7 


18 


1.5 


17 


1.1 


3 


4.5 


36 


4.4 


38 


— 



TYPICAL SEQUENCE, 


PERFECT CODE 


• 






SNR (2kHz) 


SMOOTHED OUTPUT 


VITERBI 


OUTPUT 


IMPROVEMENT RATIO 
IN BIT ERROR RATE 


NO PRE-FILTER 


BIT 


LTR 


BIT 


LTR 




(db) 


(%) 


(%) 


(%) 


(%) 




6 


0.70 


8 


0.37 


3 


1.9 


5 


1.1 


12 


0.80 


10 


1.4 


4 


1.2 


15 


0.90 


13 


1.3 


3 


2.9 


24 


2.7 


22 


1.1 
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TABLE XI 



COMPARISON OF ML AND MAP ESTIMATES 



SNR (2kHz) 

NO PRE-FILTER 
(db) 


MAXIMUM LIKELIHOOD 
ESTIMATOR 

BIT ERROR RATE 
(%) 


MAXIMUM A-PCSTERIORI 
ESTIMATOR 

BIT ERROR RATE 
(%) 


IMPROVEMENT 

RATIO 


6 


0.45 


0.37 


1.2 


5 


0.91 


0.80 


1.1 


4 


1.1 


0.90 


1.2 


3 


2.8 


2.7 


1.0 



Typical sequence 25 wpm, perfect code 
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bit error rates of about 3% or higher remain practically 
unaffected. Again, it was noted that most errors are 
isolated dot insertions during word-spaces. Since the third 
order model accounts for the small transition probability 
of .such events, such a model should virtually eliminate this 
type of error, resulting in a very significant decrease in 
error rates. 

Table XII shows a typical decoded sequence as output by 
the algorithm, with the filtered and smoothed outputs shown 
for comparison. 'This sequence was part of the 35 wpm AR 
sequence at 5 dB SNR in 2 kHz. Since it was noted that the 
algorithm created some new errors as well as correcting 
errors, likelihood computations could probably be improved 
by employing better character- length density estimation 
procedures . 

Two examples of likelihood computations and the evolution 
of the most likely sequence as the algorithm progresses are 
shown in Tables XHIa and b. Table XHIa is an example of 
decoding with highly probable likelihoods, while Table XHIb 
shows an example of a correction by the algorithm. The array 
of numbers headed by SURVIVOR SEQUENCES indicate to which 
previous node the node at stage k is connected? the length 
of each survivor sequence is shown immediately below this 
array. The line through the survivor sequence array shows 
the final path for the minimum- length sequence. The next 
array shows the computed log- likelihoods for each character. 
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followed by the minimum length sequence with its length. 

(The number 200 was sufficiently large to be used for 
infinity.) Finally, the mark and space durations as 
determined by the thresholded output are shown for comparison. 

In these tables (Xllla and b) , the following numbers 
correspond to the decoded morse characters: 

1 - dot 

2 - element-space 

3,4 - letter-space 
5,6 - dash 

The letter-space and dash are represnted by two numbers since 
such decisions may result either from comparisons made on a 
straightforward character versus character basis or from 
comparisons made utilizing the modified likelihood computation 
previously described. In each case the sequence /n./n— ^.~ 
was transmitted; Table XHIb shows how the Viterbi algorithm 
corrected the thresholded sequence ( A ~ . -n . ~ . ~) to form a 
dash from the inner . /s . sequence. This situation is similar 
to that shown previously in Figure 23. 
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TABLE XII 



COMPARISON OF VITERBI OUTPUT WITH 
SMOOTHED AND FILTERED OUTPUTS 



FILTERED 


AR 


E 


AR 


AR 


SMOOTHED 


AR 


E 


MR 


AR 


VITERBI 


AR 


E 


AR 


AR 



FILTERED 


RR 


AR 


E AR 


SMOOTHED 


AR 


AR 


AR 


VITERBI 


AR 


AR 


AR 



UC 


AL 


AR 


AR 


E 


AR 


E 


UR 


AR 


AR 


AR 


AR 


AR 


AR 


E 


AR 


E 


AR 


AR 


AR 


MR 


AR 


AR 


AR 


E 


AN 




AR 


AR 


AR 


AR 


AR 


AR 


AR 


I 


AG 










AR 


AR 


AR 


AR 


I 


AG 










AR 


AR 


AR 


AR 




AR 











Note : 

1. Sequence of AR at 35 wpm, 5 db SNR, tranmitted. 



SUMMARY OF VITERBI OUTPUT: 



New errors made 2 
Errors uncorrected 2 
Corrections 4 
Net improvement 2 
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TABLE XHIa. VITERBI DECODER COMPUTATION 
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TABLE Xlllb. VITERBI DECODER COMPUTATION 
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VII. PRE-DETECTION MODEL AND FILTER 



Since it is assumed that the receiver has been tuned to 
the signal carrier frequency, it is sufficient to model the 
signal as one of known frequency and unknown amplitude and 
phase of the form A sin(cot+0). The received signal then is 
of the form 

z(t) = A,sin(cot+0) + v(t) 

where v(t) represents additive white gaussian noise. A state 
model of the signal process may be obtained by rewriting the 
transmitted signal in the form a sin cot + b cos cot [4] . 
Letting represent the discretized signal in this form, the 
following state model results: 



x^k+l) 




1 


COT 




x 1 (k) 


x 2 (k+l) 




-COT 


1 




x 2 (k) 



where x 2 (t) “ b cos cot - a sin cot 

k = time index 
T = sampling interval 



100 



The observation model then is 



z(k) = [ 1 0 ] 



x 1 (k) 
x 2 (k) 



+ v(k) 



This model may be written more compactly in vector 
notation : 



x(k+l) = $(k) x (k) 
z (k) = H (k) x (k) + v (k) 

where is the state transition matrix and H is the measurement 
vector. The on-off keying nature of the signal may be 
accounted for in an intuitive way simply by multiplying the 
observed value, z(k), by the probability that the signal is 
present. This probability is readily obtained from the 
demodulated output of the filter by using the algorithm 
previously derived for mark and space transition probabilities 
for the baseband signal. 

B. FILTER ALGORITHM 

The general filter algorithm presented in Appendix A, 
using the above signal model, was used to filter the (down- 
converted) IF signal. Demodulation was accomplished digitally 
by squaring and averaging the x-^ state estimate. The 
subroutine used previously for the calculation of Q was then 
used in exactly the same way as for the baseband model except 
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that the indices k and j were advanced by an amount of 
time equal to the delay due to averaging in the demodulation 
process. Additionally it was found that using zero proba- 
bility during space intervals was too low to allow recovery 
when the signal pulse occurred. The probability 0.5 was 
found to be sufficient and was used whenever zero probability 
would normally have been called for. 

C. IMPLEMENTATION AND RESULTS 

This filter was programmed and tested using test signals 
of -2 dB though +2 dB SNR in a 2 kHz bandwidth, and -15 dB 
through -11 dB in a 2 kHz bandwidth with 100 Hz bandpass 
filtering prior to sampling. In order to determine the 
effectiveness of modifying z(k) by the probability as des- 
cribed above, the filter was also run with the probability 
set to 1. The signal was sampled at 4000 samples per second, 
the (down- converted) carrier frequency was 1000 Hz, and the 
modulating signal was a square wave with period equivalent 
to a code speed of 25 wpm. Because of the lengthy processing 
time (1 second of data required about 30 seconds of processing 
time) , no large sample error rates as such were obtained. 
However, as can be seen from the output. Figures 24 through 
31, the filter performed well on signals of SNR -1 dB and 
above (2 kHz) and -14 dB and above with 100 Hz pre-filtering. 

Figures 24 through 27 are typical examples of the test 
run made using the straightforward filter with no modification 



of the observed signal. In these figures channel 1 is the 
input signal, channel 2 is the filtered signal, and channel 

3 is the demodulated output. Figures 24 and 25 are the 
results obtained from the signal in the 2 kHz bandwidth; 
Figures 26 and 27 show the results of 100 Hz analog bandpass 
filtering prior to sampling. 

Figures 28 through 31 are the outputs of the processor 
using the modified observation model. Channel 1 is the 
input signal; channel 2 is the input multiplied by the 
probability obtained from the transition probability estima- 
tion algorithm; channel 3 is the filtered signal, and channel 

4 is the demodulated output. The input signals were the 
same as those shown in Figures 24-27. 

Table XIV presents the results of a bit error rate 
analysis made on a sample size of approximately 100 bits. 

Since the test signal was not morse code, no letter error 
rates were obtained. The projected letter error rates 
shown in the table were determined simply by multiplying 
the bit error rate by 10 since this is approximately the 
proportionality factor between these two error rates , as 
can be seen from Tables III-VI. 

Based on this limited data, the conclusion may be tenta- 
tively drawn that a processing scheme employing 100 Hz 
analog filtering followed by discrete optimum linear filtering 
and detection will yield acceptable decoded error rates on 
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FIGURE 24. Unmodified Pre-Detection Filter Output, 2 dB SNR, 2 kHz BW 
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FIGURE 25. Unmodified Pre-Detection Filter Output, -1 dB SNR, 2 kHz BW 
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FIGURE 26. Unmodified Pre-Detection Filter Output, 
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FIGURE 27. Unmodified Pre-Detection Filter Output, -14 dB SNR, 100 Hz BPF 
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FIGURE 28. Modified Pre-Detection Filter Output, 2 dB SNR, 2 kHz BW 
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FIGURE 29. Modified Pre-Detection Filter Output, -1 dB SNR, 2 kHz BW 
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FIGURE 30. Modified Pre-Detection Filter Output, -11 dB SNR, 100 Hz BPF 



Channel 




CM 



fO 






111 



FIGURE 31. Modified Pre-Detection P'ilter Output, -14 dB SNR, 100 Hz BPF 



TABLE XIV 



ERROR RATES FOR PRE-DETECTION FILTERING 



(a) No 100 Hz BP Filter; 



SNR 

(2 kHz) 
(dB) 



Bit Error Rate 
(%) 



Projected Letter 
Error Rate 



(%) 



2 



0 



0 



1 

0 

-1 



0 



0 



0 0 
1 10 



-2 



2 



20 



(b) With 100 Hz BP Filter: 



SNR 

(2 kHz) 
(dB) 



Bit Error Rate 

(%) 



Projected Letter 
Error Rate 
(%) 



-11 

-12 

-13 

-14 

-15 



0 

0 

0 

1 

1 



0 

0 

0 

10 

10 



NOTE: Rates are based on small sample 

size of approximately 100 
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signals of SNR as low as -14 dB in a 2 kHz bandwidth. 
Smoothing of the state estimates at this processing stage 
would probably yield additional gain although an additional 
penalty in processing time would be incurred. 

Although this scheme meets the design objectives 
specified for the processor, it is probably unrealistic to 
assume that such processing could be implemented in real 
time on any existing machine without parallel processing. 
Coding efficiency could undoubtedly be improved but the 
inherent number of multiplications involved in the filter 
algorithm would eventually limit the processing speed. 

One possible alternative is to down-convert the received 
signal to the 0-100 Hz band, and with a 100 Hz low-pass 
filter in place, sample at 200 samples/sec. With such 
filtering, the processing described above would be feasible. 
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VIII. CONCLUSIONS AND RECOMMENDATIONS 



A. CONCLUSIONS 

The results of the various processing schemes considered 
in this investigation lead to the following definite 
conclusions . 

1. Kalman filtering of the demodulated morse signal 
provides at least a 50% reduction in letter error rate over 
unprocessed letter error rates of approximately 60% or less. 

2. Smoothing of the demodulated signal provides a 
reduction in letter error rate of approximately 30% over 
filtering alone, resulting in an overall reduction of 65% 
over the unprocessed error rate. 

3. Viterbi decoding of the morse characters, using the 
smoothed state estimate for likelihood computation and using 
a first order Markov model of the code to obtain transition 
probabilities, provides a further significant reduction in 
error rate if the input error rate is on the order of 10% 

or less. For input error rates larger than 10%, Viterbi 
decoding is of little value. 

4. In almost all cases, the processor performed better 
on sloppy code and typical sequences than on perfect code 
and AR sequences. This behavior is to be expected since the 
statistics of the input signal should match those for which 
the processor was designed in order to achieve optimum 
performance . 
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5. By using a 100 Hz bandpass filter, pre-detection 
Kalman filtering allows recovery of the code with a resulting 
error rate of about 10% or less for signals whose SNR in a 

2 kHz bandwidth is -14 dB or higher. 

6. The Pickering 230-D decoder (see Appendix C) allows 
recovery of the code with a resulting error rate of about 
10% or less for signals whose SNR in a 2 kHz bandwidth is 
approximately -6 dB or higher, also by using a 100 Hz 
bandpass filter. 

By making reasonable extrapolations to various other 
processing arrangements based on the results presented here, 
the projected performance of several alternative processors 
may be obtained with a high degree of confidence. The 
processing arrangements presented in Table XV are listed more 
or less in order of increasing complexity and power, showing 
the letter error rates to be expected from each; the projected 
figures are indicated by an asterisk. The rates presented 
assume that the processors, except for the human operators, 
are preceded by a 100 Hz bandpass filter. 

Since no measurements were made for SNR's above -6 dB , 
except for the Pickering 230-D, the projected rates in the 
-4 dB column were all determined by extrapolation of the 
measured values. The projected values for coherent demodula- 
tion (row 6) are based on the presumption that such demodula- 
tion would provide approximately 3 dB of processing gain over 
non-coherent demodulation; the subsequent values for Kalman 
filtering, smoothing, and Viterbi decoding (rows 7-9) were 
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TABLE XV 



PERFORMANCE OF ALTERNATIVE PROCESSING SCHEMES 



SNR (in 2 kHz) (dB) 



Processor 


-15 


-14 


-13 


-10 


-7 


-4 


1. Nan-coherent 
demodulation 


>50 


>50 


>50 


>50 


35 


10 * 


2 . ( 1 ) followed by 
Kalman Filter 


>50 


>50 


>50 


>50 


15 


5* 


3. (1) followed by 
smoothing 


>50 


>50 


>50 


>50 


9 


3* 


4. (3) followed by 
Viterbi Decoder 


>50 


>50 


>50 


>50 


3 


1 * 


5. Pickering 230-D 


>50 


>50 


>50 


~50 


15 


7 


6 . Coherent 
Demodulation 


>50* 


>50* 


>50* 


35* 


10 * 


5* 


7. ( 6 ) followed by 
Kalman Filter 


>50* 


>50* 


>50* 


15* 


5* 


3* 


8 . ( 6 ) followed by 
smoothing 


>50* 


>50* 


>50* 


10 * 


3* 


2 * 


9. ( 8 ) followed by 
Viterbi Decoder 


>50* 


>50* 


>50* 


4* 


2 * 


1 * 


10. Pre-detection 
Kalman Filter 


20 * 


10 


5* 


0 * 


0 * 


0 * 


1L ( 10 ) followed by 
Viterbi Decoder 


15* 


4* 


2 * 


0 * 


0 * 


0 * 


12. Typical Operator 


? 


7 


7 


7 


10 * 


5* 


13. Good Operator 


7 


7 


10 


5 


1 


0 * 



Letter error rates (%) *Projected 
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obtained by applying the respective improvement ratios 
actually obtained for non-coherent demodulation. In the 
case of the pre-detection Kalman filter (row 10) , a 10% error 
rate was actually obtained for the -15 dB entry, and no errors 
were obtained for the -13 dB entry. However, since the small 
sample size tends to decrease the confidence in these figures, 
these rates were increased on a subjective basis by observing 
the quality of the output signal. Again the subsequent rates 
for Viterbi decoding (row 11) were obtained by applying the 
previously determined improvement ratio. The figures for 
the human operators are the author's estimate based on the 
error rates obtained previously for the amateur radio 
operators. A degradation factor of 3-6 dB should probably 
be added to the SNR's shown for typical field-operation 
performance . 

In summary it may be concluded that the performance of 
a good operator can be approached by a processing scheme 
employing pre-detection 100 Hz bandpass filtering followed 
by discrete Kalman filtering and a Viterbi decoder. The 
performance of a typical operator can be obtained by either 
the 230-D or by coherent demodulation preceded by 100 Hz 
bandpass filtering, with further improvement provided by 
Kalman filtering, smoothing and Viterbi decoding. 

B . RE COMMEN DAT I ON S 

A more appropriate model for the demodulated baseband 
signal would be one which incorporates the exponential 
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rise-times and fall-times of the signal pulses instead of 
the abrupt rise-times and fall-times as used here. A filter 
incorporating such a model would probably have the effect 
of reducing the errors caused by extremely short dot 
insertions during letter-spaces and word-spaces. Better 
estimators for received demodulated signal noise power and 
for character duration mean values and/or probability 
densities should be tested, although it is felt that little 
advantage would be gained from more sophisticated techniques 
in an actual operational environment. A third (or higher) 
order Markov model of the code would undoubtedly show a 
significant improvement over the first-order model in the 
effectiveness of the Viterbi decoder at sufficiently low 
input error rates. Such Markov modeling, along with improved 
likelihood calculations, deserves further investigation. 

The pre-detection Kalman filter shows the greatest 
advantage, and more investigation of such filtering, and 
possibly smoothing, with the goal of making the processor 
real-time, should be undertaken. Down-conversion of the IF 
signal to 100 Hz or so followed by a 100 Hz bandpass filter, 
and sampling at 200-400 samples/sec, may be the solution. 
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APPENDIX A 



THE DISCRETE KALMAN FILTER 

The discrete Kalman filter is an extension of the more 
familiar Wiener, minimum mean-square-error, matched filter 
to nonconstant coefficient multivariable systems with 
nonstationary noise, implemented in sequential, or recursive, 
form. The result of every processing cycle is the current 
estimate of the state under consideration. As each new 
observation is made, the current estimate is updated to 
reflect the information content of this new measurement. 

A brief outline of the derivation of the filter algorithm 
is presented below. 

The message model is described by the linear vector 
difference equation: 

x(k+l) = _£(k)x(k) + £(k)w(k) 

where the input noise, or random forcing function, w is a 
zero-mean white-noise process, with covariance 

cov[w(k) ,w ( j) ] =Q(k)S k (k-j) 

The observation, or measurement, model is given by the linear 
algebraic equation 
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{ 



:z(k) = H(k)x(k) + v(k) 

where v(k) is also a zero-mean white-noise process, with 
cov [v(k) ,v( j ) ] =R(k)6 k (k-j) . 

For simplicity, it is assumed that w and v are uncorrelated, 
and that w and the initial value of x are uncorrelated. 

The first step in obtaining the estimate at time k is 
to predict ahead from the estimate obtained at time k-1. 

This prediction may be expressed as the conditional 
expectation : 



x (k | k-1) = E [x(k) | z (k-1) ] . 



But since the measurement z_(k-l) is embedded in the previous 
estimate, x(k-l|k-l), this prediction may be expressed as 



x (k | k-1) = ^x(k-l|k-l) + _TE [w (k-1) ] . 



Since E[w(k)] = 0, this reduces to 



x (k | k-1) = ^x(k-l | k-1) . 

The estimate at time k may then be determined by 
considering it to be a summation of this prediction plus 
a correction term employing the measured value: 
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I 

where z_(k|k-l) is determined as follows: 

£(k | k-1) = E[z(k) |z(k-l) ] 

= E [ (Hx (k) + v(k) ) | z^(k-l) ] 

= E[Hx(k) |z(k-l) ] + E[v(k) |z(k-l) ] . 

Under the given assumptions, this reduces to: 

£(k | k-1) = H E[x(k) |z(k-l) ] 

= H ± x (k-1 1 k-1) . 

The filter equation (A-l) then becomes: 

x(k|k) = (k-1 1 k-1) + G(k) [z_(k) - H$x (k-1 1 k-1) ] 

where G(k) for optimum filtering is yet to be determined. 

It is desired to determine G(k) such that the resulting 



estimates are optimum in a minimum mean-square-error sense. 
Since the error in the estimate is given by: 



the covariance matrix of estimation error is 



V(k|k) = E [ (x(k) - x(k |k) ) (x(k) -x(k|k)) T ] 

Minimizing a scalar quadratic form based on this matrix 
will yield the optimum gain G(k) . Derivations are available 
in several texts [4], [11], and will not be presented here. 

The optimum gain results in the following algorithm for the 
discrete Kalman filter. 

Gain: 

G (k) = V(k|k-1) H T [H V(k|k-1)H T + R]" 1 



Estimation : 

x(k|k) = x(k|k-l) + G(k) [£(k) - Hx(k [k-1) ] 

Estimation variance: 

V(k | k) = [I - G (k) H ] V(k | k— 1) 



Prediction variance: 

V(k+l|k) = 4> V(k | k) + Q (k) 



Prediction : 



APPENDIX B 



THE VITERBI ALGORITHM 



An instructive and informative description of the 
Viterbi algorithm is presented in [8]; the following outline 
is based on this article. The MAP sequence estimation 
problem previously stated in Section VI. A. is formally the 
same as the problem of finding the shortest path through a 
certain graph, called a trellis diagram. In this diagram, 
illustrated in Figure 31, each node corresponds to a distinct 
state of the Markov signal process and each branch represents 
a transition to some new state at the next instant of time. 

In this representation, every possible state sequence 
corresponds to a unique path through the trellis. 

Each branch is assigned the length 

Me k ) = -In p ( x k+ 1 l x k ) - ln p ( z k U k ) 



where 



represents a transition from x k to x^. + ^ 



The p(x k+ii x k ) 

while P (z k | C k ) 



are the sequence transition probabilities 
are the likelihoods of a particular state. 
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The algorithm for determining the shortest path through 
the trellis diagram is a version of forward dynamic 
programming and is stated as follows: 

M = number of states; K = length of input sequence 

1. STORAGE: 
k 

x(x k ) , 1 £ x k <_ M 
£(x k ) , 1 1 x k <_ M 

2. INITIALIZATION: 
k = 0 

/\ A 

x(x^)=x q , x(m) arbitrary, m ^ x q 
r(x Q )=0; £(m) = °°, m ^ x q 

3. RECURSION: 

a. Compute: 

— ( x k+l' x k) ~ I.( x k ) ( C k ) , for a11 ^k ; 

b. Find: 

£(x k+ i) = min r( x k+1 ' x k ) / for ea ch x k+1 ? 



(time index) 

(survivor terminating in x k ) 
(survivor length) 
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c . 



Store nx k+1 ) 

A 

x(x k+ i) ; 



and the corresponding survivor sequence 



d. Set k to k+1 and repeat until k = K. 



As an example [8] , the recursive determination of the 
shortest path through the trellis shown in Figure 32 is 
shown in Figure 33. 
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k=0 



k=2 



k=3 



k=4 



k=5 



k=1 




Figure 32. Trellis Diagram with Assigned Lengths 
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Figure 33. Example of Viterbi Algorithm 
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APPENDIX C 



PICKERING 230-D PERFORMANCE EVALUATION 



The Pickering 230-D is a currently available automatic 
decoder and transcriber which includes some front-end 
processing. The front-end consists of what is essentially 
a type of coherent detection scheme with a bandwidth of 
about 180 Hz; the center frequency of the unit tested was 
at 875 Hz. The letter error rate versus SNR of the 230-D 
was determined in the laboratory using the model KB-1 
keyboard to send perfect code. The error rates were 
determined for both a 2 kHz bandwidth and a 100 Hz bandwidth 
and are shown in Table XVI. A performance evaluation of the 
230-D for error rates obtained using actual operator- 
transmitted code with no noise is available in Ref. [2] . 
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TABLE XVI 



Pickering 230-D Error Rates 



SNR 


35 


wpm 




30 wpn 


25 


wpm 


(in 2 kHz) 


2 kHz 


100 Hz 


2 kHz 100 Hz 


2 kHz 


100 


(dB) 


(%) 


(%) 


(%) 


(%) 


(%) 


(%: 


1.6 


1 


3 


1 


0 


0 


0 


0.4 


4 


15 


4 


6 


0 


0 


-0.9 


29 


6 


22 


13 


21 


6 


-2.0 


- 


7 


- 


11 


- 


6 


-4.4 


- 


7 


- 


7 


- 


4 


-6.4 


- 


15 


- 


15 


- 


5 


-8.9 


- 


21 


- 


27 


- 


22 


-10.5 


- 


60 


- 


78 


- 


61 


-12.4 


_ 




_ 


_ 


_ 


— 



Notes : 

1) Sample size: Approximately 100 letters in each case 

2) Perfect Code 
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COMPUTER PROGRAMS 



sF0RTRAM LS#G0 



1 

2 

3 

4 

5 

6 

*7 

/ 

8 

Q 

1C 

11 

12 

13 

1* 

15 

16 

17 

18 
15? 
20 
21 
22 

23 

24 

25 

26 

27 

28 

29 

30 

31 

32 

33 

34 

35 

36 

37 
3« 

39 

40 

41 

42 

43 

44 

45 

46 

47 
4 8 



cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 



c c 

C THIS PROGRAM I^FlE^ENTS THE POST- DETECTION' KALMAN C 

C FILTER A NO L I NEAR SMOOTHER/ ALONG WITH THE REOjIRED C 

C SIGNAL PARAMETER ESTIMATOR ALGORITHMS. THE following c 

C UTILITY SUBROUTINES ARE USED: C 

c buffer i n - reads data on input tape. c 

C DAL - D/A CONVERSION AND 0'JT D UT ROUTINE C 

C THE FOLLOWING SUBROUTINES ARE PROVIDED: C 

C BACK - BACKWARD FILTER ROUTINE C 

C STATS- ESTIMATES SIGNAL PARAMETERS C 

C V A R W - COMPUTES THE VARIANCE OF W ( K ) C 

C C 

C INPUT PARAMETERS AT TIME OF EXECUTION ARE I C 

C INITIAL ESTIMATE OP T C 

C THE PARAMETERS T1 AND T2 C 

C THE DOT-DURATION DENSITY PARAMETER# LDOT C 

C THE DASH-DURATION DENSITY PARAMETER/ DDASH C 

C C 

CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 

c 

REAL LDOT#MEAMX#MEAN#MEANXB# v NX1/M\|X2, m ER] t 



DIMENSION XHS M N1 ( 5 0 ) # X H S M ( 150)/XHS M N2( 1CC)#XHSMD( IOC) 
DIMENSION EVARF2(250) #XHATF2(250) #EVARF(250) 

DIMENSION X"P( 1 25 ) / PV AR?2 ( 125)#XHAT32( 125) 

D I MENS I ON I BUF ( 500 ) , XHATF ( 250 ) 

D I MENS ION XHATB1 ( 250 ) # P V ARBI ( 250 ) # X3 ( 25C ) / XHT ( 1C ) 
DIMENSION MERIT (4/ 2# 12 )# TRANS ( 1 1 ) 

CPMM0NYBL0CK1/X\T/XNTOXHSM/XH$ m D/ IU/ I J 1 / IU2/ IUL/ I UK # 
ITS# TCS/FMEANS/NSHQRT 
C0MM0N/3L 0CK2/NElE m # MFR I T 
C9M^QN/BL0CK3/ TRANS 
C 

NAMELIST LDOT# DDASH# T# Tl# T2 

c ENTER TRANSITION PROBABILITIES FOR VITPR5I DECODER 
C 

TRANS( 1 ) = - ALPG ( 0 ♦ 6 82 ) 

TRANS ( 2 ) = - ALOG ( 0 * 3 1 8 ) 

TRANS(3)=-AL0G(0«5O4 ) 

TRANS ( 4 ) = - ALOG ( 0.396 ) 

TRANS ( 5) *-ALOG( 0.25) 

TRAN5(6)=”AL 00(0.25) 

TRANS(7)=-AL0G(0.25) 

TRA\S(8)=-ALOG(0.25) 

TRAN5(9) =- ALOG (C. 684 ) 

TRANS( 10) = " A L 0 G ( 0*316) 

TRANS ( 11 ) = - AL W G ( 0 * 1 ) 

39 CONTINUE 



49 

50 

51 



C 

C 

C 



INPUT PARAMETERS 
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52 

53 

54 

55 

56 

57 

58 

59 

6C 

61 

62 

63 

64 

65 

66 

67 

6? 

69 

70 

71 

72 

73 

74 

75 

76 

77 

78 

79 

80 

81 

82 

83 

84 

85 

86 

87 

88 

89 

90 

91 

92 

93 

94 

95 

96 

97 

98 

99 

100 

101 

102 

103 

104 



OUTPUT (101) 'ENTER T , Th AMD T2 (MSEC)' 
OUTPUT (101) 'ENTER LOOT AND DDASH (MSEC)' 
I NPUT ( 101 ) 

C 

C INITIALIZE variables and indices 
c 

T = T/2 ♦ ; LDOT =LD0T/2 • ; DDASH = DDASH/2* 

T1 5 T 1/2 » ; T2 = T2/2 « 

A=LDOT/T 
B*DDASH/ ( 3. *T ) 

X3UTL = 0. ; NEAiMX 3 = 0 • b ; MEANX=Q*5 ; XT = T 

XTD=3.*T 

IX = 0 

KS1=KS2=0 
VA1=VA2=1 • 

SI =52= .5 

SS 1 = SS2 = * 4 

MNX 1 =MNX2=0 * 5 

IN1=IN2=IS1=IS2=0 

JT c UT? s JTD=JTD2=10 

JTC? = JTC= 10 

NELEM=NSH9RT=0 

IJ=I J1 = IJ2=I 1 = 1 I 1 = 1 12 = 0 

S m N 1 * S M N 2 = 0 • 

XHSM) = XHSN2 = C • 

XHSMN1(1)=0» ; XHSNM2 ( 1 ) =0 • 

I X N = C 

PV AR = T V AR* V ARX = 0 ♦ 5 

K I N0 = 0 ; K I MDE = 0 ; UIND = 0 ; JINDF=0 
G= 1 » ; XMATP”0« ; I X H L * 0 

F M I N = • 0 5 ; XT1=XT2 = T 
XTD1 = XTD2 = 3 • T 
XTC1-XTC2=T 
XTA=XTC=T 
DO ? <1=1/4 
. D9 2 K2= 1 / 2 
DO 2 K3=l/ 12 
2 M.ER T T ( K 1 / K2 / K3 ) = 0 « 

38 CONTINUE 
C 

C INPUT DATA FROM tape 

c. 

112 CALL BUFFERING 1, 1, If’UF, 500, IERR) 

1 IF ( I ERR. EG. 1) 33 Tft 1 

GO T9 (1/111,310*112), I ERR 
111 D2 999 I NDE = 1,2 

I N 0 E 1 = ( INDE-l ) *250+1 
INDE2= INDE+250 
C 

C BEGIN PR9CESSING 
C 

DO 999 INDEX= INDE1/ INDE2 
XI = FLO AT ( I BUT ( INDEX) )/2**23 



I T = T 
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105 

106 

107 

103 

109 

110 

111 

112 

113 

114 

115 

116 

117 

118 

119 

120 

121 

122 

123 

124 

125 

126 

127 

12 ? 

129 

130 

131 

132 

133 

134 

135 

136 

137 

138 

139 

140 

141 

142 

143 

144 

145 

146 

147 

148 

149 

150 

151 

152 

153 

154 

155 

156 

157 



UUU 



C 

c UNPROCESSED OUTPUT/ X9UT1 
C 

X9UT 1 =0 • 

I F ( X 1 • GE • MEANX ) X0UT1=1. 

C 

C IX IS THE RUNNING TIME INDEX 

c 

IX= ix+l 

IEUX.GT.250) IX = 1 
C 

c ixn inhibits smoothing algorithm until sufficient 

C DATA is OBTAINED (500 SAMPLES) 

c 

I X N = I X N + 1 

IF( IXN.GT.500) JXN=501 
C 

C STORE INPUT SAMPLES IN REVERSE ORDER FOR BACKWARD 
C FILTERING 
C 

I F ( I X * GT . 1 25 ) XB 1 = XBB ( I X - 1 25 ) ; GO TO 9 
XB 1 3 XB ( IX) 

X B B ( 126- IX) =XP (251- IX) 

9 X3<251-IX)»X1 

CALL STATS 
VARXF = VARX _ 

C 

C FORWARD FILTER A L 3 8 R I T m 
C 

399 GAINF=PVAR/(PVAP+VARXF) 

EVAR= ( 1 . -GA INF ) * P V A R 

XHAT=XHATP+GAINF*(X1-MEANX-XHATP) 

THRESH8LDED FILTERED OUTPUT/ X0UT2 

X9UT2 S 0 . • 

IF (XHAT‘GE.0. ) X0UT2U. 

I X H A T = - 1 

IF ( XHAT »GE*0. ) I XHAT = 1 

KIND=KINDF ; JI\D C JI\DF ; IXHL=IXHLF 
T = X T ; TD = XTD ; LD9T=A*XT ; DDASH=B*XTD 
XH=XHAT 

350 CALL VARw 

QF = Q* (FMEAN**2) 

PVAR=EVAR+OF 

K I NDF = K I ND ; JI\DF = JIND ; IXHLFMXHL 
XHATP=XP 
C 

C STORE FORWARD-FILTERED ESTIMATES AND VARIANCES FOR 
c smoothing 
C 

EVARFl=EVARF2 ( IX) j EVARF2( IX) =EVARF< IX) 

F.VARF (]X) =EVAR 
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170 

171 
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174 
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200 
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XHATF1 =VHATF2( IX) ; XHATF2 ( IX) = XHATF ( IX) 

XHATFf IX) *XHAT 
351 call back 
CALL V A R W 
09 = 0* ( F M E A N 3 ** 2 ) 

K I ND3=K 1 N D ; JI\'DB = JIND ; IXHLB=IXHL 
C 

C REVERSE STORED VALUES FROM BACKWARD FILTER FOR 
C SMOOTHING 
C 

IF( IX*GT. 125) XHTB1=XHATB2( IX-125) : GO TO 452 
XHTF?1 = XHATB1 ( I X ) 

XHATB2 ( 126- IX) =XHATE1 (251- IX) 

492 CONTINUE 

XHATBl (251- IX) =XHATBP 

IF ( IX. GT. 125) PVB1 = PVARB2( IX-125) ; GO TO 493 
PVB1 =PVAR31 (IX) 

P V A R 3 2 ( 126- IX ) =PVARB1 (251- IX) 

493 CONTINUE 

PVARB1 (251-IX) =°VARB 
XHATBP-XP 
PVARB C EVAR3+QD 
C 

C SMOOTHING ALGORITHM 
C 

WK = EVARF 1 / ( 1 . +EVARF1/PVB1 ) 

EVARS S ( ( 1 • - WK/PV31 )**2)*EVARF1 + (WK**2)/PVB1 
WRK = ( X H T D 1 / P V B 1 ) * E V A R S 
I F ( I X • ED • 230 ) W5K=0. 

XHATS=XHATEl/( 1 .+EVARF1/PVB1 ) +WBK 
IF ( I XN «LE • 500 ) XHATS = 0. ; TS = XT ; TDS = XTD 
C 

C AVERAGE AND STORE SMOOTHED ESTIMATE AT EACH POINT FOR 
C USE LATER IN LIKELIHOOD COMPUTATIONS 
C 

I I T = C*75*TS 
111=111+1 

I F ( I I1.3T.150) 111 = 1 

11=11+1 

I E ( I I • GT • 50 ) 11 = 1 

IF(II.GT.IIT) I 1 = 1 
XHS 1 = XHSMNl (II) 

XHS V M ( I I ) =XHATS 

XHSM1 = XHSM,1 + ( XHATS-XHS1 ) / ( 0»75*TS ) 

X'-+SM( III) = X H S M 1 
I I T P = 0 • 67*TDS 
112=112+1 

IF( I I2.GT.100) 112=1 

I F ( I I 2 . GT • I I T2 ) 112=1 

113=113+1 

I F ( I I 3 * GT • 100 ) 113=1 
XHS2 = XHSMN2 ( I I 2 ) 

XHSM.N?( I 12) =XHATS 
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233 

234 

235 
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241 
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XMSM2=XHSM2+ ( XHATS-VHS2 ) /( 0.67+TDS) 

XHSMD ( I I 3 ) =XHS^2 
I J“ I 11-10 

IF(JJ.IE.O) I J= 1 50+ I j 
I J 1 = I I1-5-I IT 
IF ( I Jl .LF.'O ) I Jl 3 150+ 1 Jl 
I J2= I 13-5 

I F ( I J2 • UE • 0 ) I J2 = 1 00+ I J2 
XSUM=0. 

DO 920 1*2/10 
XSU M = XS'JM + XHT ( I ) 

920 XM T ( i-l )=XHT( I ) 

XSU^sXSU-t + XHATS 
XMT ( 10) aXHATS 
XSU M =XSU9/10. 

X9UT = 0 » 

IF(XSUM.GT*0. 0000001 ) XFMjT = 1 • 

C 

c estimate t and td 

c 

C XNT IS THE MARK DURATION COUNT 
C XNTC IS THE SPACE DURATION COUNT 
C 

KND = C 

IF( (X0UT-X8UTL) .LE.-0*99) KND = 1 ; GO TO 900 
IF ( X3UT • GE * 0.99) NT = NT + 1 } GO TO 990 
NT 3 0 

GO TO 990 
900 X N T = N T 

NT =0 

IF ( (XNT»LT • 0 . 5*T 1 ) *0R» (XNT *GT * 6 » *T2 ) ) GO TO 990 
I F ( XNT * LE • ( 3 * *T 1 + T2)/2.) G9 TO 991 
JTD= JTD+ 1 

I F ( JTD.EG.10) JTD2=0 ; XTD?=XTQ1 
IFIUTD.EQ.20) JTD=1 
XTD1=XTD1+(XNT-XTD1 )/UTD 
JTD2*UTD2+1 

IF ( JTD2.EQ. 10) JTD 1 =0 ; XTD1 3 XTD2 
IF (JTD2.EQ.20) JTD2 = 1 
XTD2=XTD2+ ( XNT -XTD2 ) / JTD2 
XTD-XTDl 

IF ( JTD2.GT. JTD1 ) XTD=XTD2 

IF ( XTD.UT • ( 2. *XT ) ) *TD = ?«*XT; XTD 1 = XTD2 = XTD 
1 JTD 3 JTD?=5 
GO TO 990 
991 JT 3 JT+ 1 

IF(JT.EQ.IO) JT2 = 0 ; XT2 = XT 1 
I F ( JT . EQ ♦ 20 ) JT 3 1 
XT1»XT1 + (XNT-XT1 )/JT 
JT2 = JT2+ 1 

I F ( JT ? • EQ » 1 0 ) JT = 0 ; XT 1 =XT2 
IF ( JT2.EQ.20) JT2 3 1 
XT2*XT2+ (XNT"XT2)/JT2 
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X T A = XT 1 

IF ( JT2.GT . JT ) XT A = X T 2 
990 CONTINUE 

IF ( ( X0UT-X9UTL ) »GE«C«99) <ND = 1 i GO T^ 1000 
IF(XSUT«LE«0.01 > NTC=NTC+ 1 ; GO T9 1090 
NTC = 0 ; GO T9 1090 
1000 XNTONTC 
N T C c 0 

I F < (XNTC.LT*0*5«T1 ) .ef?» (XNTC.GT.6. *T2) ) GO T 9 1090 
IF ( XNTC ♦ LE • ( 3 • * T 1 + T2 ) /2 • ) JTC = JTC+1J G9 T9 1091 
G9 TO 1090 

1091 IF < JTC.EG.10) JTC2 = 0 ; XTC2 = XTCl 
I F < JTC.EQ.20) JTC = 1 

XTC1 = XTC1 + (X,\TC-XTC1 >/JTC 
JTC2= JTC2+1 

IF < JTC2*EG. 10) J T C = 0 ; XTC1 = XTC? 

IF ( J T C 2 • E 0 * 2 0 > JTC2 = 1 
XTC2=XTC2+ (XNTC"XTC2)/JTC2 
XTC=XTC1 

I F ( JTC2.GT. JTC) XTC = XTC2 

c 

c SENSE SWITCH 1 IS USED TO SELECT VITFRbI DECODER 
c 0 p T I 0 N DURING PROGRAM EXECUTION 

C 

1090 I F { SENSE SWITCH 1) 1092,1093 

1092 I r (NTC.GT*h.*TS) KND= 1 i XNTC-NTC i .NTC = 0 
IF(KND.Efi.l) CALL LIKELIHOOD 

1C93 CONTINUE 

XT= (XTA + XTO/2. 

X9UTL=X9UT 

C 

c output variables tq d/a routine for analog recording 

c 

TN0RN=XT/1OO. 

call DAL ( xl / XC’JTI/ X0UT2+ XHAT/XHT81/XHATS, X0JT, VARX, 
lGF/GAlNF# TN9RM) 

999 CONTINUE 
GO TO 38 

310 OUTPUT ( 101 ) 'END OF RUN ; HIT * TO GO ' 

INPUT < 101 ) 

GO TO 39 
C 
C 
C 

subroutine back 

c 

C THIS ROUTINE IS THE BACKWARD FILTER 

C 

C INITIALIZE INDICES FOR BACKWARD FILTER BASED ON 
C FINAL VALUES OF FORWARD-FILTERED ESTIMATES 
C 

IF ( I X • \'E • 1 ) GO TO 490 
IF(KINDF.EQ.C) <INDB=0 ; GO TO 360 
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K I NDB=XT -K I MDF 

IF ( K 1 ND3 « LT • 0 ) KIN9B*XTD-<I\DF 
I F ( K I \'D0 • l_T • 0 ) <INDr = 0 

360 IF< JINOF.CQ'O) UINDB=0 ; GO TO 361 
JIN0B=XT« JINOF 

IF (JIN03.LT *C) JIN35=XTD-J1\DF 
I F ( J I NDQ • LT • 0 ) U I N 0 P = 0 

361 CONTINUE 

XHAT3P = XrlAT ; I XML 5 = I XH A T 

PVAPB=EVAR j VaFXB=VaFX ; MFANXBsMEANX ; GB=GF 
T S = T 3 ; TDS*TDB ; FMEANS = FME AN'B 

fmeanb*fmean ; T3=xt ; tdb=xtd 

A 90 CONTINUE 

c 

C FILTER ALGORITHM 
C 

GA I NB=PVARB/ ( PVAR3+VARX3 ) 

EV ARB = ( 1 .-GAINS) *PVARB 
XHATB*XHATBP+GAlNB* ( X31-MEANX5-XHAT5P ) 

IXHATSa-1 

IFIXHATB.GE.O. ) I X H A T 3 = 1 
I X M A T = IXHAT3 

A91 !<INO = KINDB ; J I ND = J I NDB ) IXHL = IXHLB 

I XHLL = I XHLL3 ; T = T3 ; TD = T[5B ; LOOT = A * T J DDASH=3*TD3 
XH=XHAT5 
A50 RETURN 
C 
C 
C 

SUBROUTINE STATS 
C 

C THIS SUBROUTINE CALCULATES THE INPUT SIGNAL MEAN/ 

C VARIANCE/ AND ESTIMATES THE SIGNAL AMPLITUDE AND SNR 

C 

C 

C INPUT MEAN AND VARIANCE CALCULATIONS 
C 

I M 1 = I M 1 + 1 

IF( IM1.EQ.500) I m 2 = 0 ; VA2 = VAl i MNX2 = MNX1 

IF ( I Ml .EG* 1000) I M 1 = 1 

MNX 1 =MNX 1 + ( X 1 -MNX 1 ) / I M 1 

VA 1 =VA1 + ( < X 1 - M E A N X ) **2-VAl )/IMl 

IM2* r'2 + 1 

IF( IM2.FQ.500) I M 1 = C ; VA1*VA2 ; MNX1 = MNX2 

I F ( IM2.FQ.1000) I M2 = 1 

MNX2*mnX2+(X1-nNX2)/IM2 

VA2-VA2+ ( (Xl-MEANX) **2 -Va2)/IM2 

ME A\'X = MNX 1 

V A * V A 1 

IF ( IM2.GT • I M 1 ) MEANX=MNX2 ; VA=VA2 
C 

c signal amplitude estimator 
c 
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388 
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393 
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IF ( XI .LE.ME4NX) GO TO 80 
IS1-IS1+1 

I F ( I S 1 • F.O • 500 ) I S2 = C l S2 = S1 
IF( IS1.FR.100C) I S 1 = 1 
Si 'SI +( XI -Si I/ISl 
132=1 S2+1 

IF ( IS2.EC.50C ) IS1-0 i S1=S2 
I F ( IS?* FQ. 1000) IS2*1 

S2=S2+ (X1-S2)/1S2 
S ' S l 

IF( IS2.GT.IS1 ) S = S2 
30 TO 81 

80 CONTINUE 
KS 1 =KS 1 + 1 

IFUS1 .FQ.500) KS2 = 0 ; SS? = SSl 
IF(KS1.E0«1C00) KS1=1 
SS1=SS1+(X1-SS1 )/KSl 
<S2=KS2+1 

IF (KS2.EQ.500) KS1=0 ; SSl=SS2 
IF (XS2.EQ. 1000) KS2= 1 
,SS2=SS2+(X1-SS2) /KS2 
SS=SS1 

IF ('<S2.GT.KS1 ) SS = SS2 

81 CONTI NUF 
F‘MEAN = S-3S 

IF (FMEAN.LT • F M I \ ) F N E A N = F N I N 
C 

C NOISE VARIANCE ESTIMATOR 
C 

VARX=VA-0«9*(F M EAN**2)/4. 

I F ( VARX • UE * 0 . ) V A RX = 0*0000001 
RETURN 
C 
C 

c 

c 

SUBROUTINE V A R W 
C 

c q estimation algorithm 

C KIND'TIME INDEX/ </ for mark probabilities 

C JIND'TI^E INDEX/ J/ FOR S D ACE PROBABILITIES 

C 

C (FOR USE WITH PRE-D FILTER/ INITIALIZE KlND/JlND TO 50) 
C 

xp=xh; p«i.o 

IF ( I X H L • E Q » 0 ) GO TO 10 

I F ( ( ( IXHAT-IXHL) «GT.O) .AND. ( jind.st.t-ldot) ) JI \D=0; 
1KIN0=0J go to 10 

IF( ( IXHAT-IXUL)»GT.0) IXHAT=-IXHAT ; XP=-XH 
I F ( ( ( IXHAT-IXHL) .LT.O) .and. (KIND.GT.T~L DOT) ) J!\D=0; 
KIND = 0; GO TO 1C 

IF ( ( IXRAT-I XHL ) .LT.O) IXHAT=-IXHAT J XP=-XH 
10 IF ( IXHAT*EO. 1 ) K I \'D = K I ND+ 1 
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1FUXHAT.EG.-1 ) JIND=JIND+1 i GO T9 20 
C 

c MARK PROBABILITIES 
C 

XK = K I N-'D 

IF ( XK • LEI • (T-LD9T) ) GO T 5 11 
GO TO 1? 
li- G = 0 » 

IF(KIND»EQ»1> X P = F M E A N / 2 ♦ 

GO TO 100 

12 IF( (XK.LE. (T+LDOT) ) . AND • (XK.GT. (T-LD0T) ) ) GO TO 13 
GO T0 14 

13 P W K 0 1 = (T-XK 1/(4. *LDOT 1+0.75 
PWKN 1=1* - PwKOl 

Q = P W K 0 1 * P W K N 1 ; PsPWKOi; GO TO 100 

14 IF( (XK.LE. (TD-DDASH) ) .AND* (XK.GE. (T4LD0T) ) ) 0 = 0.; 

1G? TO 100 

IF( (XK.LE* (TD+DDASH) ) .AND* (XK*GE. (TD-DDASH) ) ) GO TO 17 
GO TO 18 

17 PWK01=(TD-XK)/(2.*DDASH) + 0.5 
PWKMsl. -PWKC1 

q=pwkoi*°wkni ; p=pwkci; go to 100 

18 G=0»25j GO TO 100 

c 

c SPACE PROBABILITIES 
C 



20 


X J = J I NO ; P = 0 • 5 

IF (XJ.LF. ( T-LD9T 1 1 00 TO 21 






GO TO 22 




21 


0 = 0. 

I F ( J I NO * EG • 1 1 XP = -FMEAN/2. 
GO TO 100 




22 


IF ( ( XJ.LE* (T + LDOT 1 1 .AND* ( X J . GE . ( T -LOOT 1 1 1 
GO TO 24 


GO TO 23 


23 


PWK00=(-XJ/(2»*LDBT)+0.5*( l.+T/LDOT) 1*12* 
PWKNO= 1 . -PWKOO 

g=pwkoo*pwkno; p=pwkoo+o.5; GO TO ICC 


/17.+5./17 


24 


I F ( (XJ.LE. (TD-DDASH) 1 .AND* (XJ*GE. (T+L00T1 1 1 3*0.; 
1 GO TO 100 



IF( (XJ.LE* (TD+DOASH) ) .A V D* (XJ.GE. (TO-OOASH) 1 ) GO TO 27 
GO TO 28 

27 PWKOO 3 ( ( -XJ/( 2. « DO ASH) +0.5* ( 1 • +TD/DDASH 1 11*0.8+0.2 
PWKNO= 1 • -PWKOO 

q=pwkoo*pwkno; p=pwkoo+o.5; go to 100 

28 IF (XJ.LE.5.U 1 0 = 0. ; GO TO ICO 
PWKOO =25. «T**2/XJ**2 

P W K N 0 = 1 .-PWKOO; Q = PaK00*PWK\0; P = PWK00 + 0 * 5 
ICO CONTINUE 

IF(XK«GT.10.*T) KIN[5 = 0 ; 0=C.?5 
IF (XJ.GT . 10* *T ) J I NJD = 0 ; 0 = 0.25 
IF (P.GT. 1 . 1 P=1 .0 
1XHL= IXMAT 

RETURN 
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COMPILATION! 

SUBROUTINE LIKELIHOOD 



C 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

r 

c 



THIS subroutine CALCULATES THE LIKELIHOODS OF EACH 
CHARACTER , USING THE M ETHoD GIVEN IN THE TEXT OF 
THE THESIS 



MERIT { I# liN) 
MERIT (2/ 1/N) 
MERIT ( 3 > 1 > N ) 



MERIT ( 3/ 2/ N ) 



MERIT (4. 1 # N ) 
MERIT ( 4 / 2 » N ) 



IS THE DOT LIKELIHOOD 
IS THE ELEMENT -S°ACE LIKELIHOOD 
IS THE LETTER-SPACE LIKELIHOOD EOR 
COMPUTATIONS MADE ON A CHARACTER VS 
CHARACTER BASIS 

IS THE LETTER-SPACE LIKELIHOOD E9R 
COMPUTATIONS INVOLVING THREE SHORT 
CHARACTERS 

IS THE CHAR VS CHAR DASH LIKELIHOOD 
IS THE THREE SHORT CHAR DASH LIKELIHOOD 



REAL MERIT 

DIMENSION MERIT! 4>2> 12)/XMERT (3) 

DIMENSION -XHSM( 150); XHS m D< 100 ) aESTORE ( 8 ) 

DIMENSION MARK ( 11 )*NSPACE< 11 ) 

C9Mv<0N/6L0CK1/XNT>X\TOXHSMjXHSMD/ I J, I J 1 a IJ2/ IJL* UK/ 
ITS/ TD5,F v EANS>NSH0RT 

c 9 M m o N / B L e C K 2 / N E l E m / M E R I T 

C0MM0N/BL0CK4/^ARK# NSPACE/ NCHAR 

NCHAR=NCHAR+1 

MARK (NCHAR) =2 . *X.NT 

NSPACE ( NCHAR ) =2* *XNTC 



C 

C MAKE TENTATIVE DECISIONS BASED ON THRESHOLDED 
C S M 69 THED OUTPUT 

C 

IF ( X N T . L T • 0 • 5 ) GO TO 103 

IE (XNT •GE»2»*TS) EDuR = XNT / NSH0RT=O; ES'JM = 0.; 

1 GO TO 102 

I F ( (XNT«LT.2.*TS) .AND* (XNT»GT»0»5) ) EDUR = XNT; 

1G9 TO 101 

103 JF(XNTC«GF.2.*TS) edur=xntc; NSH0RT=O) esum*o»; 

1 GO TO 102 

I E ( ( XN'TC.LT *2 *0*TS) - AND - ( XNTC • GT . 0 . 5 ) ) EDUR = XNTC 



C 

C COMPUTATIONS FOR SHORT CHARACTERS 
C 

101 NELEM = NELEM + 1 

NSHORT = NSH9RT + 1 

IF ( NSH9RT *GT . 8 5 NSHQRT=0 
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81 
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EST0RF<NSH0RT>*EDUR 
ESU'-sESUM + EDlJR 

I F ( M5H0RT » GT . 3 ) ESJN=£SUM-EST0RE ( NSHORT-3 ) 

TM1* 1 . 

IF (FDURfGT.TS) TNI = 2. - ( EDUR + TS ) / ( 2, +TS ) 

TM2 = 0 • 

IF ( EDUR • GT • TDS/3 • ) TM2-1 ( TDS-EDUR ) / ( 0 . 67 *TDS ) 

X M EF I Ts ( XHSM ( I J) /FNE A\'S + 0. 5 ) *TM1 * ( 1 . - T M 2 ) 
(XVCRIT.LT.C. ) xnlkit=o. 

IF ( X M ERI T .GT. 1 . ) XMERIT=1. 

MERIT ( 1 , 1 ,NCl_EM) = X V ER I T 

XMER I T = ( 1 . - ( XHS^ ( I J ) /FME ANS+O • 5 ) ) * T" 1 * ( 1 * - TM2 ) 

I F ( XN-ER I T • L T • 0 • ) X^ERIT=0. 

IFtXNERlT.GT.l.) XMERIT=1. 

MER I T ( 2i IM.ElCM) =Xmt RI 7 

,X M ER I T = ( X H 5 M ( I J)/FMEA\S + C.5)« ( 1 . -TM1* ( 1 . -TM2 ) ) 

IF (XNERIT.LT.O. ) X m £ R I T = 0 • 

IF (XMERI r.GT.l. ) X M E R I T = 1 • 

M ER I T ( 4 / 1 1 NElEM ) aXNTRlT 

XNER I T = ( 1 . - ( XHSM ( I J ) /FMEAN5 + Q » 5 ) > « ( 1 . - TMlM 1 . - TM2 ) ) 
IFIXNERIT.LT.O. ) XMFRIT=0. 

IF ( XMERI T .GT . 1 . ) XMERIT=1. 

MERIT(3/ 1, NELEM) »XMERIT 
C 

C MODIFIED LIKELIHOOD COMPUTATIONS FOR 3 OR MORE 
C S u ORT CHARACTERS IN SUCCESSION 
C 

I F ( NSHORT » L T ♦ 3 ) GO TO 300 
I F ( NELE M * L T • 3 ) GO TO 300 
DO 200 M/3 

200 XMERTI I )*AMAX(MERIT(1/1jNELEM-3 + I ), 

1MERIT(2> 1/NELEV-3+I ) ) 

X N R T = ( l . -XMERT ( 1 ) ) * ( 1 . -XkERT ( 2 ) ) * ( 1 • ~X V ERT ( 3) ) 
XMRT=XMRT**0. 33333 
TM2= l . 

IF(ESIJM.lT.TDS) TM2=l.-(TDS-ESUM)/(0,67+TDS) 

IF (TM2.LT. 0. ) T M 2 * 0 . 

X V EP I T = ( XHSMQ ( I J2)/FMEANS+0*5 >*TM2*XMRT 
IF (XMERIT .LT.O. ) XMFRIT*0. 

IF(XVERIT.GT.l.) XMER 1 T= 1 • 

M E R I T ( 4 1 2/ NELEM-2 ) = XMER I T 

XMERIT* ( 1 »-(XHS M .D( I J2 ) /FME ANS + O » 5 ) )«T m 2+XMRT 
I r ( XMER I T . LT . 0 . ) XMERIT=0. 

I c ( XMER I T « GT . 1 . ) XMERIT*1. 

MERI T(3>2/ NELEM-2 )= XMERIT 
GO TO 300 
C 

C COMPUTATIONS for long CHARACTERS 

c 

102 NELEM-NELCM+3 
I JK= I J-EDUR/3. 

I F < I JK . LE • 0 ) I JK = 1 50+ I JK 
IJL=IJ”(2.*EDUR)/3. 
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101 

102 

103 

104 

105 

106 

107 

108. 

109 

110 

111 

112 

113 

114 

115 

116 

117 

118 

119 

120 

121 

12 ? 

123 

124 

125 

126 

127 

128 

129 

130 

131 

132 

133 

134 

135 

136 

137 

138 

139 

140 

141 

142 

143 

144 

145 

146 

147 

148 

149 

150 

151 

152 

153 



IF ( I JL.L.E .0) IJt = 150+IJL 
T^1=0* 

I F ( CUUR • l_T * 3 . * TS ) m=2.-(EDUR + TS)/(2.*TS) 

I F ( T M 1 • L T • 0 « ) T M 1 * 0 • ; I F ( T M 1 • GT • 1 ♦ ) T M 3 = 1 • 

TM2 = 1 . 

I F ( FOUR • L T • TRS ) TM2=l.-( TDS-EDUR ) / ( C . 67 « TDS ) 

X M E R 1 T = ( XHSMD ( I J2) /FM,EANS+0»5 ) *TM2* ( 1 . « TM1 ) 

IF (X.-TRIT.LT.O. ) XMER'IT = 0. 

I F ( X v F. R I T . G T « 1 ♦ ) XMER I T = 1 • 

MERIT! 4/2/ NELEM-2>= XMER IT 

XMER I T = ( 1 . - ( yHS M D ( I J2 ) /FhEANS+ 0.5 ) ) *TM2« ( 1 #-TMl > 

IF < XYER I 7 .LT.C. ) XMER1T = 0. 

IF(XFERIT'GT.l. ) XMERIT=1» 

MERI T( 3/2/NEL pM -2) = XMER I T 
C 

C MODIFIED COMPUTATIONS FOR DOT AND ELEMENT-SPACE 
C LIVELIHOODS 
C 

IF(E0UH/3»*LT»2.5*TS) 30 TO 3C0 

XMER I T s'( XMSM( I JL ) /FMEANS + O * 5 ) * ( 1 . -TM2* ( 1 • *TM1 ) } 

IF (X V FR!T • L T . 0 . ) XMER1T = 0* 

IF ( XMERI T .GT . 1 . ) XMER1T=1* 

MERIT ( 1/ 1jNElC m -2)sXMER!T 

XMEP I T= ( 1 . - ( X^S"' ( I JL) /F M EANS + C.5) ) * ( 1 . - TM2* ( 1 . - T M 1 ) ) 
IF ( XMERIT .LT.C. ) XMFR I T - C • 

IF(X“ERIT.GT.1. ) XMERIT»1* 

M E R I T ( 2 , 1 , N E L E M - 2 ) = X m e R I T 

XH ER I T = ( XHSM( I J< ) /FvLANS + 0.5) * ( 1 . -TM 2* ( 1 • -T^ 1 ) ) 

I F ( XMER I T *LT . 0 . ) XMERIT = 0. 

IF(XMERIT.GT.l. ) X M F R I T = 1 • 

MERIT (1,1 , NELEM-1 )= XMER If 

X m E R I T = ( 1 ♦ - ( XHSM ( I J< ) / F y EANS + 0»5) ) * ( 1 . -TM2* ( 1 « - T M. 1 ) ) 
IF(XMERIT.LT,0. ) XMERIT=0. 

I F ( XMER IT . GT • 1 • ) X M E R 1 T = 1 * 

MERIT (2, 1, NELEM-1 )= XMER IT 

XMER I T= ( XHSM ( I J)/F^EA\S + 0.5)*( 1 ♦ -TM2* ( 1 . -TMl ) ) 
IF(XMERIT*LT.O. ) X y ERIT=0* 

IF(X-'ERIT.GT. 1 . ) X m C R ! T = 1 1 
MERIT! 1/ 1 >.NELEM) =XMErit 

X y ER I T = ( 1 • - ( XH$ y ( I J ) /F.MEANS + O • 5 ) ) * ( 1 . -T^'2* ( 1 . - TMl ) ) 

1 r ( XMER I T • LT t 0 • ) X M E R I T = 0 • 

IF ( XMER I T .GT. 1 . ) X M ERIT=1« 

MERIT (2, 1 / NELE M ) = X M E R I T 
300 XNT C = XNT = 0* 

IF ( NELEM. GE *9) GO To 301 
RETURN 
C 
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154 

155 

156 

157 

158 

159 

160 

161 

16? 

163 

164 

165 

166 

167 

163 

169 

170 

171 

172 

i 



ir<M ERIT( I 1 > i 3, I 2 ) » LE • 0 • ) MERIT < I 1/ 13/ 12) =200* ; 
1 GO TO 501 

M ER I 7 < 11/ 13 / l2)=-AL0G( k 1FRIT( II, 13 / I?) ) 

5C1 CONTI fun 

CALL DECODER 
C 

C INITIALIZE LIKELIHOODS TOR NEXT ITERATION 
C 

nchar=o 

DO 600 1=1/11 

m ar'< ( i ) =njsrace n ) =0 

600 CONTI NUF 

DO 505 11=1/4 
DO 505 12=1/12 
DO 505 13=1/2 
505 HPRIT ( I 1/ 13/ 12) =0. 

NELEM=0 

RETURN 

END 
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9F 

1 

2 

3 

4 

5 

6 

7 

R 

9 

10 

11 

1? 

13 

14 

15 

16 

17 

18 

19 

20 

21 

22 

23 

24 

25 

26 

27 

28 

29 

30 

31 

32 

33 

34 

35 

36 

37 

38 

39 

40 

41 

42 

43 

44 

45 

46 

47 



COMPILATION 

SUBROUTINE DECODER 
C 

C THIS SUBROUTINE IMPLEMENTS the 
C V1TER81 ALGORITHM 
C 
C 

REAL LSAV/MERI T»LAMOA /LENGTH, NU 

REAL LNGTH1 /LNGTH2/ L N G T H 3 / LNGTH4, LNGTH5, LNGTH6 
INTEGER SEQSAV/CHAR 

DIMENSION LSAV(6/12), ISAVE(6/12)/WSAVE(5) 
DIMENSION TRANS( 11 ) /XTRANSI 11 ) 

D I MENS 1 ON LENGTH ( 6 ) / Y L N 0 T H 16 ) t WLNG7 h ( 6 ) / 2LNGTH ( 6 ) 
DIMENSION MER1 T ( 4 / 2 / 1 ? ) / CH AR ( 6/ 1 2 ) 

DIMENSION MAQK( 11)/ NSPACE ( 11 ) 

DIMENSION NSEG( 12)/ NSQ< 12) /SEQSAVC5/ 12 ) 
C0MMSN/BL0CK2/NELEM/ MER I T 
C9MM0N/PL0C<3/ TRANS 

co mm on/block 4 /mark/ nspace/ nchar 

DATA L/0/ 

c 

C INITIALIZATION 

c 

DO 30 1=1/12 
80 NSEG ( I ) =NSQ ( I ) =0 
C 

C USE A PRIORI PROBABILITIES OF EACH CHARACTER 

C FOR INITIALIZATION 

C 

LENGTH ( 1 ) =-ALOG(0.269) 

LEiNGTH(B)=-AL?G(0.341) 

LENGTH ( 3 ) = LENGTH ( 4 ) = - AL R G ( 0 . 1 59 ) 

LENGTH ( 5 ) =LENGTH ( 6 ) =~ ALOG ( 0 *232 ) 

DO 100 1=1/6 

NLNGTH ( I ) =ZLNGTH( I ) =Yl.NGTH( I ) =200. 

ICO CONTINUE 

DO 101 1=1/11 
XTRANSI I ) =0. 

101 CONTINUE 

DO 10 1=1/6 
DO 10 J= 1 / 1 2 
LSAV ( 1/ J)=0. 

I SAVE ( 1/ J)=0 
10 CONTINUE 

c 

C SENSE SWITCH 2 IS USED TO SELECT OUTPUT OPTION 

C 

IF (SENSE SWITCH 2) 8/9 
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48 

49 

50 

51 

5? 

53 

54 

55 

56 

57 

53 

59 

60 

61 

62 

63 

6 4 

65 

66 

67 

68 

69 

70 

71 

7? 

73 

74 

75 

76 

77 

73 

79 

80 

81 

8? 

83 

84 

85 

86 

87 

88 

89 

90 

91 

92 

93 

94 

95 

96 

97 

98 

99 

100 



8 WRITE (6*3000) 

3000 F 9 R M A T ( ' 1 ' ) 

WRITE (6/ 4010) 

4010 F0R M AT(5X* 'SURVIV9R SEQUENCES’) 
WRPE(6*4000) (1/1 = 1/ 11) 

4000 FORMAT* 14Xz 12( I2/5X) ) 

9 CONTINUE 
C 

c MAIN ALGORITHM 

c 

C N n DE 1 IS A DOT 
C NODE 2 IS AN ELEMENT SPACE 
C NODES 3 AND 4 ARE LET TER-SF'ACES 

C NODES 5 AND 6 ARE DASHES 

C 
C 

NELM1=NFLEM+1 
DO 900 < = 1 / UELM 1 
IF ( K • EG • 2 ) GO TO 910 
GO TO 911 

910 DO Pll J=l/ll 
XTRANS(J) = T R A N S ( J ) 

911 CONTINUE 

IF U.HE.NELNl ) GO TO 920 
DO 920 J=l/ll 
XTRANS(vJ) =0. 

920 CONTINUE 
C 

c NODE 1 

LAMQA=LENGTH(2)+XTRANS(3) ; 1=2 

NU a LbNGTH(3)*XTRANS(5) 

IF (NU*LT. LAMDA) lamca = nlj ; 1=3 

NU = ZL-NGTH ( 4 ) +XTRANS ( 5 ) 

IF(NU*LT. LAMDA) LAM[)a = NU ; 1=4 

I SAVE ( 1/0 = 1 

LN3TH1=LAMDA+MERIT ( 1/ 1/0 

300 CONTINUE 
C 

C NODE 2 

LAMDA = LENoTH( 1 ) +XTRANS* 1 ) ; 1 = 1 

NU=ZLNGTH( 4 ) +XTRANS ( 6 ) 

IF (NU«LT. LAMDA ) LAMDA=NU i 1=4 
NU=LENSTH (5 ) +XTRANS ( 9 ) 
IF(nU‘LT.la m da) lamda=nu ; 1=5 
NU'ZLNGTH ( 6 ) +XTR ANS ( 9 ) 

IF(NU«LT. LAMDA) LAMDA=NU } 1=6 

I SAVE ( 2/K ) = I 

LNGTH2 = L AMDA + ^ER I T (2 / 1/0 

301 CONTINUE 
C 

C NODE 3 

LAMDA=LENGTH( 1)+XTRANS(2) ; 1=1 
NU=LENGTH(5)+XTRANS( 10) 
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1C1 

1C? 

103 

104 

105 

1 06 

107 

108 

109 

110 
1 1 1 
11 ? 
113 
1 1 4 

115 

116 
117 
113 

119 

120 
121 
122 

123 

124 

125 

126 

127 

128 

129 

130 

131 

132 

133 

134 

135 

136 

137 

138 

139 

140 

141 

142 

143 

144 

145 

146 

147 

148 

149 

150 

151 

152 

153 



I E ( NU • LT • LA'iO A ) LAMDA=NU ; 1=5 

N'J = ZLNGTH<6)+XT r ?AN3( 10) 
IF<\'U*LT»LAMDA ) LAMDA = NU ; 1=6 

NU»LENGTH(2)+XTRANS( 11 ) 
IF(NU»LT.LAMDA) LA v !DA = N'J ; 1=2 

I SAVE (3/0 = 1 

LNGTh3=LANDA+MERIT(3/ 1/0 

302 CONTINUE 
C 

c NODE 4 

LANDA =LENGTH < 1 ) +XTR ANS ( 2 ) i 1 = 1 
\U = 7.LNGTM ( 4 5 +X TRAMS ( 7 ) 

I F ( \U ♦ l_T • L A M OA ) LAMOA = \'U J 1=4 
NU=LEN3TM(5 )+XT«AMS (10) 

I F ( NU ♦ LT • LA MU A ) |_^*DA*NU i 1*5 
NU = ZLN'GTH ( 6 ) + XTR ANS ( 10) 

I F ( NU • LT • LAND A ) LAM[)A = NU ; 1=6 

NU = LENGTH(2)+XTRA\'S( 11 ) 

IF ( \U*LT. LAUDA ) l.A^DA = N'J ; 1=2 

ISAVE(4/K)=I 

LNGTH4=LAMDA+MERIT( 3/2/K) 

303 CONTINUE 
C 

C NODE 5 

LANDA = LENGTH( 2)+XTRANS( 4 ) i 1=2 
NU = LENGTH ( 3 > +XTRANS (8) 

IF ( NU • LT • LAND A ) LAMDA=NU i 1=3 
NU s ZLNGTM(4 )+XTRANS<8) 

I F ( NU • LT « L AMpA ) LAMDA = NU ; 1 = 4 

I SAVE (5/ 0 = 1 

LNGTH5=LANDA+^ERIT ( 4/ 1/ K ) 

304 CONTINUE 
C 

C NODE 6 

LANDA=LENGTH(?) +XTRANS( 4 ) ) 1=2 

NU = LENGTH (3 > +X TRANS (8) 

IE ( NU»LT • LAMOA ) L A N 0 A = N U i 1=3 
NU=ZLNGTH ( 4 ) +XTRANS ( 8 ) 
IE(NU*LT.LAMDA) LAMDA=NU ; 1=4 

I SAVE (6/0 = 1 

LNGTH6=LA V DA+MER I T ( 4/ Z, K ) 

305 CONTINUE 
C 

C STORE LENGTHS FOR SURVIVOR SEQUENCES 
C 

LENGTH ( 1 ) =LNGTHl 
LENGTH ( ? ) =LNCjTH2 
LENGTH! 3) =LNGTH3 
LENGTH ( 4 ) = L N G T H 4 
LENGTH (5)=LNGTH5 
LENGTH ( 6 ) =LNQTH6 
DO 11 1=1/6 

11 LSAV( I/O =LENGTH( I ) 
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15/* 

155 

156 

157 

155 

15? 

160 

161 

162 

163 

16 z* 

165 

166 

167 

165 

169 

170 

171 

172 

173 

174 

175 

176 

177 

175 

179 

180 

1S1 

182 

183 

184 

185 

186 

187 

188 

189 

190 

191 

192 

193 

194 

195 

196 

197 

198 

199 

200 

201 

202 

203 

204 

205 

206 



DO 200 1=1/6 
ZLNGTH ( I ) = Y l M G T H ( I ) 

Y L N 0 T H ( 1 ) sWLNGTU ( I ) 

WLNGTH ( I ) 'LENGTH ( I ) 

200 CONTINUE 
9C0 CONTINUE 
C 

C DETERMINE MINIMUM LENGTH SEQUENCE FROM THE 
C SIX SURVIVOR SEQUENCES 
C 

NU=LENGTM( 1 ) i ISURV'l 
D9 901 1=2/6 

I E ( LENGTH ( I ) .LT • NU ) NU = LENGTH ( I ) ; IS'JRVM 

901 CONTINUE 

XEIGHT=LENGTH( ISURV) 

I M = 1 

y,s o 

N= I SAVE ( IS'JRV/NELNl ) 

NSEC ( 1 ) =N 
32 CONTINUE 

I r ( ( N • C Q • 4 ) • 0 R • ( M ♦ E G • 6 ) ) 36 TO 30 

M=M+1 

IF ( (NElHI-M) .LF.*0) GO TO 35 
NS = I SAVE ( N/ NELMl - M ) 

GO T 6 31 

30 M = M + 3 

I E ( ( N E L N 1 - M ) • L E * 0 ) GO TO 35 
N5= I SAVE <N/ NELMl -M) 

31 I M = 1 M + 1 

NSEG ( IM) = N = N S 
GO T 8 32 
35 CONTINUE 

DO 40 1 = 1/ IM 
40 NSG ( I ) =NSEQ < IM-I+1 ) 

L n L + 1 

WSAVF.(L)=wEIGHT 
DO 71 1 = 1/ IM 
71 SEQSAV ( L/ I ) =NSG ( I ) 

I M 1 = IM+1 
DO 73 I = I M 1 / 1 2 
73 SEQSAV ( L/ I 5=0 
C 

C OUTPUT 
C. 

IF (SENSE SWITCH 2) 60/70 
C 

C OUTPUT OPTION 1 - DETAILED OUTPUT 
C 

60 CONTINUE 

DO 12 J= 1 / 6 

12 WRITE (6/ 4C01 ) J/ (ISAVE( J/0/K = l/NELCM) 

4001 FORMAT (SX/ 12( I 1/6X> ) 

WRITE (6/ 4002) 
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207 

208 

209 

210 

211 

21 ? 

213 

214 

215 

216 

217 

218 

219 

220 

221 

222 

223 

224 

225 

226 

227 

228 

229 

230 

231 

232 

233 

234 

235 

236 

237 

23S 

239 

240 

241 

24? 

243 

244 

245 

246 

247 

24 S 

249 

250 

251 

25? 

253 



4002 FORMAT ( / ) 

WR I TE ( 6/ 5000 ) 

5000 FORMAT (5X/ 'SURVIVOR LENGTHS') 

09 13 J= 1 / 6 

13 WRITE (6/4003) J/ < LSAV ( J/ K ) / <= 1 / NELEM ) 

4003 FORMAT (8X/ I1/PX/12(F6»2/1X)) 
wR I TE ( 6/ 4002 ) 

WR I TE ( 6/ 1 003 ) 

1003 E9R M AT(5X/ 'LIKELIHOODS' ) 

IJ=1 
JI =0 

D9 903 1=1/2 
JI = JI+1 

903 /.RITE (6/ 1002) JI / (MERI T ( 1 / I J/ K ) / K*1 / NELEM ) 

DO 905 1=3/4 

,D9 9C5 I J = 1 / 2 
JI sJl+1 

905 WRITE (6/ 1002) J I / ( MER I T ( 1 / 1 J/ < ) / K = 1 / NELEM ) 

1002 F6RMAT( 8X/ I 1/ 2X/ 12 (F6.2/ IX > ) 

WRITE (6/4002) 

WRITE (6/ 4021) 

4021 FORMAT ( 63X/ 'DEC90ED* ) 
wR I TE ( 6/ 4020 5 

4020 F3R m AT ( ?4X/ ' M ARK ' / 3X/ ' S U ACE ' / 12Xz 'LENGTH' /9X, 

1 • SEQUENCE ' > 

904 WRITE(6# 1000) MARK ( 1 ) / NSPACE ( 1 ) / WE I GMT/ < NSQ ( K ) / 

1 K = 2 / I M ) 

1000 FORMAT (20X/ I5/3X/ I5/11X/F9.3/ 1CX/ 1211 ) 

DO 96C I=2/NCHAR 

WRITE (6/ 1000) MARK! I )/ NSPACE ( I ) 

960 CONTINUE 

I r ( L • F CJ • 5 ) L = 0 
GO TO 72 

70 IF ( L * N E « 5 ) GO TO 72 

C 

C OUTPUT OPTION 2 - DECODED SEQUENCES AND LENGTHS ONLY 

C 

L = 0 

WRITE (6/ 6000) U'SAVE( I )/ 1 = 1/5) 

WRITE (6/ 6001 ) ( ( SEQSAV ( I/J)/J=2/12)/l=l/5) 

6000 FORMAT (5X/5(F9»3/2X)//) 

6001 FORMAT (5X/5( 11 Il/lX) ) 
wR I TE ( 6/ 6003 ) 

6003 F a R v AT(//) 

72, CONTINUE 
RETURN 
END 
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^FORTRAN LS/GO 



C 

C 

C 



THIS PROGRAM 
DESCRIBED IN 



implements 

THE THESIS 



THt PRE-DETECT ION filter 

TEXT 



REAL LDOT/MCANX/ i dent 

DIMENSION I DENT (2/2) /GAIn!2)/H( 2) 

DIMENSION TEMP M ( 2/ 2 ) / TEMPM 1 ( 2/ 2 ) / TE^PR ( ? ) / TEMPC!2) 



EVAP<2/2 



DIMENSION 
DIMENSION I3UFJA000) 
NAMELIST LDOT/DDASM/T/Tl, 



) / PVAR ( ?/ 2 ) / PH I (2/2)/ PH IT (2/ 2) 



T2/FREQ 



10 

11 

12 

13 

14 

15 

16 
17 
IS 



C 

C 

C 



OUTPUT (101) 'ENTER T AND FREQ' 
OUTPUT (101) 'LOOT AND DD ASH ' 
INPUT! 101 ) 

initialize variables AND INDICES 



XP 1 = X D 2 = 0 • 5 
TAU=0. 00025 
KT = 0 



19 

20 
21 
22 

23 

24 

25 

26 

27 

28 

29 

30 

31 

32 

33 

34 

35 

36 

37 

38 

39 

40 

41 

42 

43 

44 

45 

46 

47 
A8 

49 

50 

51 



KIND = JI\'D = 0 
XOUT=C. 

MEANX=0. 

ss = o. 

VARx=0«5 

T = 4 . * T ; LD0T*4.*L00Ti DDASH = 4 . «DDASh 

T1=4.*T1 ; T2=4 ♦ *T2 

FREQ=6.?8319«FREQ 

TD = 3 • * T 

I XHL = 0 

P=1.C 

DO 20 1=1/2 
PHI ( 1/ I )=1. 

PH I T ( I / I ) = 1 • 

20 CONTINUE 

PHI ( 1/ 2) a PH IT (2/ 1 ) =rREG*TAU 
PHI (2/ 1)=PHIT( l/2) = -FREC*TAU 
DO 21 1=1/2 
DO 21 J-l/2 
I DENT ( I / J) =0. 

EV AR ( 1/ J)*0. 

PVAR( I / J) *0* 

21 CONTINUE 

PVARM /I )=PVAR(2/2) = 1. 
IDENT(1/1)=IDENT(2/?)=1. 

H( 1 ) *1 . 

H( 2 ) =0. 

IOC CONTINUE 
C 

C INPUT DATA ERO M TAPE 
C 

112 CALL BUFFERIN'! 1/ 1/ IBUF/ 4000/ IERR) 

1 IF( IERR.EQ*!) GO TO 1 
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52 

53 

54 

55 

56 

57 

58 

59 

60 

61 

62 

63 

64 

65 

66 

67 

68 

69 

70 

71 

72 

73 

74 

75 

76 

77 

78 

79 

80 

81 

82 

83 

84 

85 

86 

87 

88 

89 

90 

91 

92 

93 

94 

95 

96 

97 

98 

99 

100 

101 

102 

103 

104 



GO TO (1/111*310/112) IERR 
111 DO 999 I\DLX=l/4000 
C 

C BEGIN PROCESSING 
C 

Xl - f-lOAT ( 1 3UF ( INDEX ) ) /2»*?3 

X2 c X 1 

X1=P*X1 

IF (XT. GT. 1000000 > KT = C 
< T = X T + 1 
CALL STATS 
C 

C GAIN 

CALL NVMULT(PVAR/H/ TEMPO 
CALL WMULT (H/ TEMPO TEMP ) 

TEM°1=1 ./(TE^p+VARX) 

CALL MVMULT (PVAR/M/ TE M PC ) 

G A I N ( 1 ) = TEMPO ( 1 ) *TEMP1 
GAIN (2 ) =TEMPC (2 ) »TEMP1 
C 

c ESTIMATION 

XH1 = XP1+GAIN( 1 ) * (X1-XP1 ) 

XH2«XP2 + GAIN(2)* (Xl-XPl ) 

C 

C ESTIMATION VARIANCE 

CALL VMULT (GAIN/ H> TEMPM) 

DO 200 1=1/2 
DO 200 J= 1 / 2 

2C0 TEMP''" ( 1/ J) = IDCNT(I/ J) -TEMPM ( I/J) 

CALL MMULT(TE v PM/PVAR/EVAR) 

C 

C PREDICTION VARIANCE 

CALL MMULT (EVAR/PHIT / T E v P M ) 

CALL MMULT(PHI /TEMP'"/PVAR) 

C 

c ONE-STEP PREDICTION 

XP1=XH1+PHI (1/2) *XH2 
XP2 = PH I (2/ 1 ) * X H 1 + X H 2 
C 

c SQUARE FILTERED ESTIMATE AND LOW-PASS FILTER 
C 

X=XH1**2 
XKT=KT • 

I F ( XT « GT « 50 ) GO TO 60 
X9UT = XOUT + ( l./X<T)*(X-xeuT) 

GO TO 61 

60 X0UT = XOUT+O«O?*(X-X0IJT) 

61 CONTINUE 
IXHAT=-1 

IF(X9UT.GE. M EANX) I XHAT = 1 
CALL V A R W 
C 

C OUTPUT VALUES TP D/A ROUTINE FOR ANALOG RECORDING 
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105 

106 

107 

108 

109 

110 

111 

112 

113 

114 

115 

116 

117 

113 

119 

120 

121 

12 ? 

123 

124 

125 

126 

127 

123 

129 

130 

131 



C 

XXHAT=FU0AT( IXMAT) 

call 0AL(X2/X1/XH1/ X0UT'GAIN< 1)/EVAR<1/1 LEVARI?/?)/ 
1VARX/P/XXHAT) 

999 CONTINUE 
310 GO TO ICO 
C 
C 
C 

subroutine stats 

c 

c t w is subroutine estimates the noise varaince 

c 

c 

XKT=KT 

IF { «T»GT .4000 ) GO TO 70 

MEANX = MFANX+ ( 1 * /XKT ) * (XOUT-MEANX) 

GO TO 71 

70 k 'EANX = MFANX+ ( 1./4C0Q. ) * ( XOUT-MEANX ) 

71 CONTINUE 

IF(X0UT*3T. m lANX) GO TO 80 
SS = SS + 0.01* (XOUT-SS ) 

80 CONTINUE 
VARX=SS 

IF ( VARX • LE • 0 « 01 ) VARX = 0.01 

RETURN 

END 
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^FORTRAN ls< go 

l: SUBROUTINE MMULT(A,e,C) 

?: c MULTIPLY TW6 matrices 

3: DIMENSION A(2/2)/8(?/2)/C<2/2) 

a: oe 10 1-1/2 

5: W 10 J--1/2 

6: ic c(I/j)=o. 

7: DC 20 1=1/2 

8: DO 20 J»li2 

9: D9 20 K*li2 

10: C( 1/ J)=A( I/K)*3(K/ J)+C( I/vJ) 

11: ?C CONTINUE 

12: RETURN 

13: END 
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9 
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COMPILATION 

SUBROUTINE MvMULT(AfBiC) 

C MULTIPLY MATRIX BY COLUMN VFCTOR 
DIMENSION A(2/2)/8(2)/C<2) 

DO 10 1=1/2 
10 C(I)=0. 

DO 20 1=1/2 
DO 20 K = 1 / 2 
C( I ) = A( I/K) *B(K)+C( I ) 

20 CONTINUE 
RETURN 
END 
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CQMPILAT ion 

SUBROUTINE VmM'JLT(A/6/C) 

C MULTIPLY ROW VECTOR 3Y matrix 
DIMENSION' A(2)/?(2/2)/C(2> 
DO 30 1*1/2 
10 C ( I ) =C » 

DO 20 1=1/2 
DO 20 K*l/2 
C( I ) = A ( K 5 * 3 ( '< / I )+C( I) 

20 CONTINUE 
RETURN 
END 
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END 0F COMPILATION 

: SUBROUTINE VV^ULTf A , E? , C ) 

: C MULTIPLY R3W VECTOR BY COLUMN VECTOR 
: DIMENSION A<2)i5(2) 

: C=0 . 

: DO 10 1=1*2 

: C c A ( I ) « B ( I ) +c 

: 10 CONTINUE 

: RETURN 

: end 
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END 0F compilation 

i: subroutine: vmult<a,B/C> 

: c multiply column vector by row vector 
: DIMENSION A(2)/3(2)*C(2/2) 

: DO 10 1=1 >2 

: DO 10 J=1;2 

: CU/J)=A( I )*B< J) 

: io continue 

; RETURN 

: end 
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